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SUMMARY 


An  earlier  effort  (reference  1)  proved  the  feasibility  of  a  procedure 
for  predicting  the  nonuniform  induced  velocity  distributions  and  the 
corresponding  airloads  experienced  by  rotor  blades.  The  present 
effort  was  undertaken  to  improve  the  aerodynamic  representations 
used  in  this  earlier  work,  improve  the  computational  efficiency,  and 
"close  the  loop"  to  enable  the  airloads  and  corresponding  blade  motions 
to  be  predicted  simultaneously. 

The  aerodynamic  representation  was  improved  by  satisfying  the 
chordwise  boundary  condition  (instead  of  using  the  lifting  line  repre¬ 
sentation).  Considerable  simplification  was  incorporated  into  the 
wake  representation  which  resulted  in  a  substantial  increase  in  the 
computing  efficiency  and  no  loss  of  accuracy  in  the  results.  The 
equations  of  motion  for  the  blade  flapping  and  flapwise  bending 
degrees  of  freedom  were  incorporated,  and  an  iterative  procedure 
was  developed  which  yields  a  simultaneous  solution  for  the  aero¬ 
dynamic  loads  and  dynamic  response  experienced  by  the  rotor 
blades. 

The  lift  loading  and  bending  moments  were  computed  for  the  H-34 
rotor  at  =0.  18  and  /£  =  0.  29,  and  the  HU-1A  rotor  system  at 

-  0.  08  and  /J.  -  0.  26.  Comparisons  of  these  results  with  measured 
results  are  presented  both  as  distributions  of  the  harmonic  components 
and  as  time  histories.  Results  of  computations  to  investigate  the 
sensitivity  of  the  computed  airloads  and  bending  moments  to  variations 
in  seme  of  the  wake  parameters  are  presented  and  discussed. 


FOREWORD 


The  work  described  in  this  report  was  performed  by  the  Cornell  Aero¬ 
nautical  Laboratory,  Inc.  (CAL),  Buffalo,  New  York,  for  the  U.  S.  Army 
Aviation  Materiel  Laboratories,  (USAAVLABS),  Fort  Eustis,  Virginia, 
over  the  15  month  period  beginning  30  March  1964.  This  work  is  an 
extension  of  work  accomplished  at  CAL  o/er  a  period  of  years  beginning 
in  I960  under  Contract  DA  44-177-TC-69B. 

Mr.  R.  A.  Piziali  was  project  engineer  and  author  of  this  report. 

Mr.  J.  Balcerak  and  Mr.  R.  Whalen  assisted  in  programming  the 
computational  procedure  for  the  IBM  7044  EDP.  The  contributions  of 
Dr.  H.  Daughaday  and  the  many  discussions  with  Messrs.  R.  P.  White, 
Jr.  ,  F.  A.  DuWaldt,  A.  R.  Trenka,  and  Dr.  I.  C.  Statler  of  CAL  were 
helpful  and  are  appreciated  by  the  author.  Mr.  J.  Yeates  administered 
the  project  for  USAA.VLABS. 


v 


BLANK  PAGE 


CONTENTS 


Page 

SUMMARY  iii 

FOREWORD  v 

LIST  OF  ILLUSTRATIONS  viii 

SYMBOLS  x 

INTRODUCTION  1 

ASSUMPTIONS  2 

WAKE  REPRESENTATION  4 

BLADE  REPRESENTATION  » 

COMPUTATIONAL  PROCEDURE  15 

COMPUTER  PROGRAM  19 

COMPUTED  RESULTS  AND  COMPARISON  WITH  MEASURED 
RESULTS  20 

CONCLUSIONS  AND  RECOMMENDATIONS  31 

REFERENCES  33 

DISTRIBUTION  70 

APPENDIXES 


I.  Two -Dimensional  Oscillating  Wing  Investigation  of  72 
Wake  Shed  Vortex  Representation 

H.  The  Chordwise  Induced  Velocity  Distribution  Due  to  81 
A  Wake  Vortex  Element 

III.  Efficient  Method  for  Numerically  Computing  gg 

Derivatives  of  Periodic  Functions 

IV.  Flow  Diagram  for  the  Computer  Program  88 

(Parts  1  and  2) 


!■  "  '  1 

1 

\ 

- - - 

.  — 

ILLUSTRATIONS 

Figure 

Page 

J 

Example  of  Possible  Wake  Configurations 

35 

2 

The  Tip-Path-Plane  Oriented  Coordinate  System  Used 
for  Computing  the  Induced  Velocities 

36 

3 

Example  of  the  /<  -Subscript  Notation  for  Airload  Points  36 

In  Rotor  Disk  When  a//?  -  3,  NA  -  8 

4 

Coordinate  System  Used  in  Blade  Representation 

37 

5 

Relation  Between  Wake  Vortex  Strengths  and  the  NRA 
Bound  Vortex  Strengths  When  tJR  =  3,  NA  -  8 

37 

6 

Major  Steps  in  the  Iterative  Procedure  of  Solving  the 
Equations 

38 

7 

Segmentation  of  Rotor  Blades  for  Computations 

39 

8 

Measured  and  Computed  Lift  and  Bending  Moments; 
HU-1A  at  /ju  =  0.  26 

40 

9 

Measured  and  Computed  Lift  and  Bending  Moments; 
HU-1  A  at  -  0.  08 

44 

10a 

Computed  Quasi-Stea.dy  F  ,  Total  F  ,  and  Lift  at 
r/R  =  0.  95;  HU-1A  at  ju,  =  0.  08 

48 

10b 

Computed  Induced  Velocity  at  r/ R  -  0.  95;  HU-1A  at 

/A  =  0.  08 

49 

11 

Measured  and  Computed  Lift  and  Lending  Moments; 
H-34  at/*  =  0,  29 

50 

12 

Measured  and  Computed  Lift  and  Bending  Moments; 
H-34  at/*  =  0.  18 

54 

13a 

Computed  Quasi-Steady  r  ,  Total  F  ,  and  Lift  at 
r/R  =  0.  95;  H-34  at  /*  =  0.  18 

58 

13b 

Computed  Induced  Velocity  at  r/R  =  0.  95;  H-34  at 
/*  =  0.  18 

59 

1 4a 

Computed  Pitching  Moments;  HU-1A  at  /c  =  0,26 

60 

14b 

Computed  Pitching  Moments;  HU-1A  at  /*  -  0,  08 

61 

viii 

Page 


Figure 


14c  Computed  Pitching  Moments;  H-34  at  ;.o  -  0.  29  62 

1 4d  Computed  Pitching  Moments;  H-34  at  /c  =  0.  18  63 

15  Computed  Pitching  Moments  about  Midchord;  H-34  64 

at  ll  -  0.  18 

16a  Comparison  of  Computed  Lift  at  r/R  -  0.95  for  65 

Wake  with  and  without  a  Concentrated  Root  Trailing 
Vortex;  HU-1A  at  /x  -  0.  26 

16b  Comparison  of  Computed  Bending  Moment  at  r/R  =  0.  21  65 
with  and  without  a  Concentrated  Root  Trailing  Vortex; 

HU-1  A  at  /x  =  0.  26 

17  Distorted  Tip  Vortices  on  Retreating  Side  of  Rotor  66 

Disk  as  Used  in  Computation;  Blades  in  the  ^  =  90° 

and  270°  Positions;  HU-1A  at  /a  =  0.  03 

18  Comparison  of  Lift  at  r/R  -  0.  95  Computed  with  67 

and  without  Wake  Distortion;  HU-1A  at  ^  =  0.08 

19  Comparison  of  Lift  at  r/R  =  0.  95  Computed  for  Two  67 

Positions  of  the  Tip  Vortex  in  the  Wake;  HU-1A  at 

/*  =  0.  08 

20  Comparison  of  the  Induced  Velocity  at  r/R  -  0.  95  for  68 

Two  Radial  Positions  of  the  Tip  Vortex  in  the  Wake; 

HU-1  A  at  jx  =  0.08 

21  Relative  Orientation  of  the  Blade  and  Recent  Tip  69 

Vortices  on  Advancing  Side  of  the  Disk;  HU-1A  at 

/x  =  0.  08 

APPENDIXES 

I,  1  Two-Dimensional  Oscillating  Airfoil  Lift  and  Pitching  77 

Moment  Transfer  Functions;  Zero  Wake  Advance  (D=  Rci) 

I.  2  Two-Dimensional  Oscillating  Airfoil  Lift  and  Pitching  7  9 

Moment  Transfer  Functions;  0.70  Wake  Advance  (Dzl-tiZ/ 3) 

IV.  1  Flow  Diagram  for  Part  I  of  Computer  Program  90 

IV.  2  Flow  Diagram  for  Part  2  of  Computer  Program  93 


1 


-M' 


SYMBOLS 


ADV 

wake  advance 

Glauert  coefficients 

AtJ  ,Sn  Fourier  coefficients  of  generalized  airload  in  the  5th  degree  of 
-  freedom 

Q-ns,b„s  Fourier  coefficients  of  blade  response  in  the  .^degree  cf  freedom 
BA  blade  azimuth  position  number 

b  blade  semichord 

airfoil  lift-curve  slope 

Cs  coefficient  of  the  quasi -steady  damping  term 

D  normalized  (by  semichord)  distance  from  airfoil  midchord  to 
first  shod  vortex 

d  normalized  (by  semichord)  shed  vortex  spacing 
fs(r)  normalized  blade  mode  shape  in  the  sth  degree  of  freedom 
Gs  generalized  airload  in  the  s** degree  of  freedom 
g  structural  damping 

fi  blade  section  plunging  velocity  (positive  upward) 

2  quasi -steady  part  of  r 

K  subscript  denoting  collocation  points  in  rotor  disk 

KA  azimuth  position  number 

£  blade  section  lift  per  unit  span 

£  stalled  blade  section  lift  per  unit  span 

&£  cross-flow  drag  force  for  the  blade  section 

Ms  generalized  mass  in  the  5th  mode 

7n  blade  section  pitching  moment  about  midchord 

ftj  stalled  blade  section  pitching  moment  about  midchord 


x 


ms(r] 

NA 

ND 

NR 

Nw 

NRA 

Pij 

b 

R 

r 

R 

Sn 

s 

t 

At 

U 

* 

v(3) 

ir(9) 
V t 

Vr 

w 


bending  moment  distribution  in  the  s*hmode 

number  of  azimuth  positions  used  in  the  computation 

number  of  degrees  of  freedom  used  in  the  computation 

number  of  blade  radial  segments  used  in  the  computation 

number  of  revolutions  of  wake  used  in  the  computation 

total  number  of  collocation  points  in  rotor  disk;  NRA  =  NR  *  NA 

wake  element  end  points 

tip  deflection  in  the  normal  modes 

total  blade  radius 

radius  to  a  blade  section 

ratio  of  computed  r  to  maximum  value  (stalled  value)  allowed 

cosine  components  of  the  chordwise  induced  velocity  coefficients 

{n  =  0,  1,  2,  .  . . ) 

subscript  denoting  each  blade  mode 
time 

time  increment  corresponding  to  each  azimuth  increment 
i.  e. ,  At  =  A  •ff/j'l 

wake  transport  velocity  in  the 2*  direction 
rotor  translational  velocity 

geometric  or  qu<.si-steady  part  of  the  chordwise  normal  velocity 
distribution 

total  chordwise  normal  velocity  distribution 

tangential  component  of  total  local  velocity  at  blade  section 

normal  component  of  total  local  velocity  at  blade  section 

total  local  velocity  at  blade  section 

wake  transport  velocity  in  the  Z-direction 


xi 


us(9)  chordwise  normal  induced  velocity  distribution 

&(/•,  ip)  specified  induced  part  of  Z- component  of  wake  transport  velocity 
X  number  of  shed  vortices  per  cycle  of  the  airfoil  motion 


chordwise  coordinate  (positive  aft) 

effective  angle  of  attack  of  blade  section  relative  to  VT 

geometric  angle  of  attack  of  blade  section  relative  to  shaft  plane 
(plane  normal  to  shaft) 

compliment  of  shaft  angle  relative  to 
tip-path-plane  angle  relative  to  In¬ 
stall  angle  for  airfoil  section 
total  bound  vorticity  of  blade  section 
chordwise  bound  vorticity  distribution 

secondary  chordwise  coordinate  defined  by  &=  -b cos  Q 

,  VrCOSCCy 

advance  ratio;  u.  =  - — 

fLr 

air  density 

induced  velocity  coefficients  of  P-  equations 
displacement  of  blade  section  relative  to  shaft  plane 
azimuth  angle 

azimuth  increment,  ^ 'J- 

NA 

blade  rotational  speed 
natural  frequency  of  the  5^mode 


Xll 


INTRODUCTION 


The  U.  S.  Army  Aviation  Materiel  Laboratories  (USAAVLABS)  hav^ 
been  directing  a  unified  experimental  and  theoretical  research  effort 
which  has  for  its  objective  the  improvement  of  the  state  of  the  art  used 
in  the  design  and  development  of  V/STOL  aircraft.  As  a  participant 
in  this  effort,  Cornell  Aeronautical  Laboratory,  Inc.  (CAL)  has  been 
working  01  the  development  of  a  more  accurate  means  (suitable  for  use 
in  the  design  and  development  of  rotor  systems)  of  predicting  the 
performance  oi  rotors  and  the  aeroelastic  response  experienced  by 
the  rotor  blades. 

The  feasibility  of  the  overall  approach  to  predicting  the  nonuniform 
inflow  distribution  resulting  from  the  wake  of  the  rotor  blades  was 
demonstrated  in  an  earlier  initial  effort  (reference  1)  conducted  at  CAL 
in  which  the  wake  configuration  was  specified  and  the  prediction  of  the 
blade  response  was  circumvented  by  assuming  the  response  to  be  given, 
This  "open  loop"  (with  respect  to  the  blade  response)  approach  enabled 
concentration  of  the  effort  on  the  aerodynamic  aspect  of  the  problem. 
Because  of  the  encouraging  results  of  the  initial  effort,  the  present 
second  phase  was  undertaken  to  improve  the  aerodynamic  representa¬ 
tion  of  the  blades  and  wake,  compute  the  aerodynamic  moments,  "close 
the  loop"  to  enable  the  blade  motions  to  be  predicted,  and  simplify  the 
computer  program. 
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ASSUMPTIONS 


The  method  developed  to  predict  the  aerodynamic  loads  and  dynamic 
response  experienced  by  the  blades  of  a  rotor  system  in  steady-state 
translational  flight  is  based  on  the  following  principal  assumptions: 

1.  The  v/ake  configuration  can  be  adequately  prescribed, 

2.  The  blade  slopes  in  the  spanwise  direction,  and  the 
section  angles  of  attack  below  stall  are  small. 

3.  The  inplare  components  of  the  induced  velocities  at  the 
tip-path  plane  are  small  and  can  be  neglected. 

4.  Below  stall,  the  lift-c.urve  slope  is  constant. 

5.  The  blade  section  circulation  is  limited  to  a  maximum 
value  for  angles  of  attack  at  and  above  stall. 

6.  For  angles  of  attack  above  stall,  the  blade  section 
normal  force  is  taken  as  the  sum  of  the  stall-limited 
circulatory  force  and  a  cross-flow  drag  force. 

7.  The  Mach  number  and  Reynolds  number  effects  are 
assumed  to  influence  only  the  lift-curve  slope. 

8.  The  interference  effects  of  the  rotor  hub,  fuselage, 
etc.  ,  are  negligible. 

The  first  assumption  (that  of  being  able  to  prescribe  adequately  the  wake 
configuration)  ic  believed  to  be  one  of  the  major  sources  of  discrepancy 
between  !;he  computed  and  rr  easured  results.  Computations  which 
indicate  the  sensitivity  of  the  predicted  results  to  the  wake  configuration 
have  been  made  and  will  be  presented  in  a  latter  section  of  this  report. 

"Wake  configuration"  means  the  time  history  of  the  spacial  distribution 
of  the  vorticity  which  streams  from  each  blade.  In  reality,  a  continuous 
sheet  of  shed  and  trailing  vorticity  streams  from  each  blade  with 
strengths  proportional  to  the  time  and  radial  rates  of  change  of  the 
bound  vorticity.  Under  the  influence  of  the  forces  of  viscosity  and  the 
three-dimensional  induced  velocity  field  (due  to  itself  and  the  blade) 
through  which  it  moves,  this  vortex  wake  diffuses,  dissipates,  and 
distorts.  The  exact  configuration  of  this  wake  is  not  completely  known 
and  the  state  of  the  art  is  not  sufficiently  advanced  to  permit  detailed 
prediction  of  this  wake  configuration.  *  The  tip  vortices,  however,  are 

*  There  is  some  (as  yet,  unpublished)  work  proceeding  at  CAL  which  has 
developed  a  method  for  computing  the  displacement  time  history  of  a 
system  of  segmented  concentrated  vortices  of  finite  core  streaming 
from  the  rotor  blade  tips  when  their  strength  distribution  is  specified. 
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known  to  exist  and  persist;  the  form  of  these  vortices  has  been  visualized 
(e.  g.  ,  the  smoke  studies  of  references  2  and  3). 

It  was  desired  to  concentrate  effort  on  the  rotor  aeroelastic  prediction 
problem  without  having  to  undertake  the  equally  complex  problem  of 
predicting  the  wake  configuration.  It  was  felt,  moreover,  that  the  wake 
configuration  could  be  prescribed  well  enough  to  allow  the  development 
of  the  remaining  aspects  of  the  method  to  proceed.  This  approach  has 
enabled  the  method  of  predicting  rotor  blade  aeroelastic  response  to  be 
developed  and  proved  practical.  The  computed  results  cannot,  of  course, 
be  any  better  than  the  prescribed  wake  configuration,  but  the  problem  of 
adequately  predicting  the  wake  configuration  can  now  be  undertaken 
separately. 

The  method  for  predicting  rotor  aeroelastic  response,  as  presently 
developed,  has  several  limitations  which  are  really  additional  assure^  - 
tions  arbitrarily  imposed  only  to  expedite  the  development  of  the  over¬ 
all  procedure;  they  are  as  follows: 

1.  Only  the  equations  of  motion  for  the  flapwise  degrees  of 
freedom  were  included.  The  equations  of  motion  for  the 
remaining  degrees  of  freedom  (i.  e.  ,  pitch,  torsion,  and 
the  inpiane  degrees  of  freedom)  and  all  the  significant 
mass,  elastic,  and  geometric  coupling  between  all 
degrees  of  freedom  remain  to  be  included. 

2.  The  rotor  control  settings  are  assumed  known,  and  the 
total  resultant  force  and  moment  computed  for  the  rotor 
were  not  required  to  balance  the  fuselage  forces  and 
moment;  i.  e.  ,  the  rotor -fuselage  trim  constraints  are 
not  imposed. 

Now  that  the  overall  procedure  has  been  developed  and  proved  feasible, 
the  additional  degrees  of  freedom  and  coupling  of  item  1,  above,  can  be 
incorporated  and  the  additional  aspect  of  the  problem  described  in  item 
2  can  be  included  to  make  the  method  complete  and  practical  for  use  in 
the  analysis,  design,  and  development  of  rotor  systems. 
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WAKE  REPRESENTATION 


As  in  reference  I,  the  shed  and  trailing  vorticity  distributions  in  the 
wake  of  each  blade  are  represented  by  an  arrangement  of  concentrated 
straight-line  vortex  segments;  however,  the  arrangement  of  these  vortex 
segments  is  now  somewhat  different  than  was  used  previously. 

In  the  present  representation,  the  azimuthal  extent,  behind  each  blade, 
of  the  grid  of  straight-line  vortex  segments  representing  the  shed  and 
trailing  vorticity  distributions  can  be  truncated  where  desired  and  the 
wake  continued  on  as  segmented  root  and/or  tip  trailing  vortices.  (The 
grid  can  still  be  carried  for  the  full  extent  of  the  wake,  if  desired. ) 

The  radial  positions  of  the  root  and/or  tip  trailing  vortices  in  the  wake 
can  be  adjusted  to  account  for  the  contraction  of  the  wake.  This  wake 
representation  can  better  represent  the  actual  wake  which,  apparently, 
rolls  up  very  quickly  into  a  tip  vortex  (see  the  smoke  studies  of  refer  ¬ 
ences  2  and  3).  Figure  1(a)  illustrates  an  example  of  a  portion  of  a  wake 
representation  with  a  root  and  tip  trailing  vortex,  and  Figure  1(b)  one 
with  only  a  tip  vortex  (for  this  example,  five  radial  blade  segments  have 
been  used). 

In  the  wake  representation  of  reference  i,  a  segmented  shed  vortex  was 
deposited  in  the  flow  at  each  azimuth  position  of  the  blade;  thus,  the 
nearest  shed  vortex  to  the  blade  was  behind  it  a  distance  proportional  to 
the  time  it  takes  the  blade  to  traverse  one  azimuthal  increment.  The 
adequacy  of  this  representation  of  the  wake  shed  vorticity  was  investi¬ 
gated  by  using  it  in  a  computational  model  to  predict  the  lift  and  pitch¬ 
ing-moment  transfer  functions  (for  both  the  pitching  and  plunging  cases) 
for  a  two-dimensional  oscillating  airfoil  at  zero  mean  angle  of  attack. 

The  results  were  compared  with  the  classical  analytical  solutions  for 
this  problem  and,  as  was  expected,  the  agreement  was  not  very  good.  A 
computational  investigation  of  discrete  shed  wake  representations  was 
then  conducted  by  use  of  this  computational  model  for  the  two-dimensional 
oscillating  wing  problem.  It  was  found  that  reasonable  agreement  could 
be  obtained  in  the  reduced  frequency  range  of  interest  by  using  the  equally 
spaced  discrete  shed  vortex  representation  if  the  entire  wake  were 
advanced,  with  respect  to  the  blade  producing  it,  a  distance  proportional 
to  70  percent  of  a  time  increment.  Thus,  on  the  basis  of  this  result,  the 
wake  representation  presently  used  is  such  that  the  wake  of  each  blade 
can  be  advanced  any  desired  amount.  A  description  of  the  two-dimension¬ 
al  oscillating  wing  investigation  of  the  discrete  shed  vortex  wake  repre¬ 
sentation  and  the  results  of  it  are  presented  in  Appendix  I. 

It  should  be  noted  tha*-  +he  two-dimensional  unsteady  wake  effects  of  thin- 
airfoil  theory  arise  f'  —i  the  shed  vcrticity;  whereas,  for  a  three-dimen¬ 
sional  oscillating  wing,  there  is  an  additional  source  of  the  unsteady  wake 
effects  which  is  the  attached  trailing -type  vorticity.  The  time -varying 
strength  of  the  trailing  vorticity  from  an  oscillating  three-dimensional 
wing  influences  the  wing  in  the  same  manner  as  the  shed  vorticity;  i.  e.  , 
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it  attenuates  and  shifts  the  phase  of  the  lift  response  to  the  time-varying 
angles  of  attack  just  as  the  shed  vorticity.  It  is  believed  that  the  repre¬ 
sentation  of  the  trailing  vorticity  reproduced  these  unsteady  wake  effects 
quite  well. 

In  reality,  the  shed  and  trailing  vorticity  distributions  in  the  wake  stream¬ 
ing  from  each  blade  are  respectively  proportional  to  time  and  radial  rates 
of  change  of  the  bound  vorticity.  Therefore,  in  the  wake  representation 
used,  the  strengths  of  the  shed  vortex  elements  deposited  in  the  wake  at 
each  azimuth  position  are  made  equal  to  the  change  in  strength  (between 
azimuth  positions)  of  the  corresponding  bound  vortex  elements.  The 
strengths  of  the  trailing  vortex  elements  (in  the  grid  portion  of  the  wake) 
are  made  equal  to  the  change  in  strength  of  adjacent  bound  vortex 
elements.  Beyond  the  truncation  of  the  wake  grid,  the  root  and/or  tip 
trailing  vortex  strengths  are  made  equal  to  the  maximum  value  of  the 
radial  distribution  of  the  bound  vorticity  on  the  blade  in  the  azimuth 
position  from  which  the  elements  were  shed.  The  effects  of  viscous 
dissipation  are  neglected  and,  thus,  the  strengths  of  the  wake  vortex 
elements  do  not  change  with  time. 

Each  v/ake  vortex  segment  end  point  is  given  an  individual  transport 
velocity  when  it  is  deposited  in  the  wake;  thus,  the  transport  velocities 
of  the  wake  element  end  points  which  originate  from  each  point  in  the 
rotor  disk  can  be  different.  All  wake  element  end  points  move  with 
their  individual  velocities  only  for  a  specified  (in  the  program  input) 
number  of  time  increments;  for  the  remaining  time  increments,  they  all 
move  with  a  common  velocity.  The  initial  wake  element  transport 
velocity  distribution  over  the  disk,  the  specified  number  of  time  incre¬ 
ments  they  move  with  their  individual  velocities,  and  the  common 
velocity  they  all  acquire  after  this  specified  number  of  time  increments 
are  all  input  parameters.  With  this  means  of  wake  description,  it  is 
possible  to  (1)  have  all  the  wake  element  end  points  move  with  a 
common  velocity  for  their  entire  life  (i.  e.  ,  an  undistorted  skewed  helical 
wake),  (2)  have  all  the  wake  element  end  points  move  with  individual 
velocities  for  their  entire  life  (i.  e.  ,  continuously  distorted  from  the 
skewed  helix),  or  (3)  have  all  the  wake  element  end  points  move  with 
individual  velocities  only  for  a  specified  length  of  time  and  then  revert 
to  a  common  velocity  (i.  e.  ,  distorted  from  the  skewed  helical  wake  for 
a  fixed  length  of  time  and  then  maintain  that  distorted  configuration 
while  moving  away  at  uniform  velocity). 

This  flexibility  of  describing  the' wake  was  built  in  only  to  allow  various- 
wake  configurations  to  be  tried  and  to  determine  the  sensitivity  of  the 
results  to  each  variation;  it  does  not  permit  true  v/ake  distortions  to 
be  used.  In  reality,  each  element  of  the  wake  experiences  a  nonperiodic 
but  time-varying  transport  velocity.  The  actual  configurations  of  the 
wake  streaming  from  each  blade  are  not  completely  understood.  For 
example,  there  is  some  question  as  to  whether  the  shed  vorticity  persists 
long  enough  or  is  in  such  a  configuration  that  it  can  significantly  influ¬ 
ence  the  following  blades.  Likewise,  does  the  trailing  vorticity  from 
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the  inboard  portion  of  the  blades  (that  vorticity  of  opposite  sense  to  that 
of  the  tip  vorticity)  persist  and  assume  a  configuration  (trailing  root 
vortex)  such  that  it  can  have  a  significant  influence  on  the  following 
blades  ? 

The  tip  vortices  and  their  configurations  have  been  made  visible  (refer¬ 
ences  2  and  3).  The  observed  distortions  of  these  tip  vortices  from 
skewed  helices  are  primarily  the  result  of  the  mutual  induced  velocities 
between  successive  coils  (layers)  of  the  wake.  These  observed  distor¬ 
tions  are  the  initial  phases  of  the  rolling -up  of  the  helices  into  two  wing- 
like  trailing  vortices  (i.  e. ,  downstream  from  the  disk,  the  rotor  wake 
is  practically  the  same  as  a  wing  wake).  These  distortions  are  influ¬ 
enced  by  the  advance  ratio,  the  disk  loading,  and  the  number  of  blades 
in  the  rotor  generating  the  vortices.  A  skewed  stack  of  equally  spaced 
parallel  vortex  rings  will  exhibit  the  same  characteristic  distortions. 

The  initial  phase  of  the  rolling -up  of  the  tip  helices  into  two  wing-like 
trailing  vortices  can  best  be  observed  in  the  downstream  views  of  the 
rotor  wake  made  visible  in  the  smoke  studies  of  reference  3. 

In  order  to  predict  the  distortions  of  these  tip  vortex  helices,  it  will  be 
necessary  to  evaluate  the  total  induced  velocity  at  each  point  on  these 
vortices  at  successive  increments  of  time  allowing  the  vortex  elements 
to  move  with  their  respective  instantaneous  velocities  during  each 
time  increment.  For  sufficiently  small  time  increments  and  enough 
collocation  points,  it  should  be  possible  to  predict  these  distortions 
(i.  e. ,  displacement  time  histories)? 

The  computation  of  the  wake-induced  velocities  at  the  blades  is  carried 
out  in  a  nonrotating  tip-path  plane  oriented  (i.  e.  ,  XY  -  plane  is  parallel 
to  tip-path  plane)  coordinate  system  which  is  illustrated  in  Figure  2. 

The  coordinate  system  translates  with  the  rotor  system  at  velocity,  Vf  ; 
thus,  the  wake  element  transport  velocities  in  the  X  and  Z  directions 
are,  respectively,  U  =  Vf  cos  ccT  and  W  =  Vf  sir,  CCr  +ur(r,y)  where  us-(r,$) 
is  a  distribution  over  the  rotor  disk  of  the  induced  part  of  the  Z  -compon¬ 
ent  of  the  wake  transport  velocity.  In  this  computational  procedure,  it 
is  this  component  which  makes  it  possible  for  each  wake  element  end 
point  to  have  an  individual  transport  velocity.  The  wake  transport  in  the 
Y  -direction  is  taken  to  be  zero;  in  reality,  the  transport  in  the  Y  -direc¬ 
tion  would  be  due  to  only  the  induced  velocities  in  that  direction  and 
would  result  in  a  lateral  distortion  of  the  wake. 

The  wake  element  end  points  in  the  wake  of  a  blade  are  labeled  as  pij 
where  L  denotes  the  azimuthal  distance  of  the  point  in  the  wake  behind 
the  blade  and  j  denotes  the  radial  wake  position;  this  notation  is  illustrated 
in  Figure  1(a).  The  equations  for  the  X  ,Y  and  Z  coordinates  of  the 
wake  points  p -y  are  (for  i  >  2) 


*  See  footnote  on  page  2. 


X  =  (Vf  cos  CcT)^~J(i-/.0  -ADV)  +  rj  cos(BA-i+A  DV)A  f  +  i1.5)bjsin  (8A  -t  +ADV)Af 


Ys  rj  sin  (BA-i+  ADV)A*  -  (1. 5) bj  cos  (BA  - i  -h  ADY) 

Z  -  -(VfSin  CCy  -  ur  (rj ,  ^ka)) ~ ~AD'/)  +■  4*00  (*j) 

where  BA  is  the  azimuth  position  of  the  blade  shedding  the  wake,  ADV 
is  the  fraction  of  a  time  increment  that  the  wake  is  advanced,  rj  is  the 
radius  to  the  blade  segment  end  points,  bj  is  the  blade  semichord  at  the 
blade  segment  end  points,  Vxa  is  the  azimuth  position  of  the  blade  when 
the  wake  element  end  point  was  shed,  and  <poe  Crj)  is  the  steady  deflection 
of  the  blade  segment  end  points.  The  coordinates  of  the  points  ptJ  when 
4  =  1  are  simply  the  coordinates  of  the  blad  trailing  edge  at  the  ends 
of  the  blade  segments. 
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BLADE  REPRESENTATION 


For  the  blade  representation,  the  theory  of  two-dimensional,  unsteady, 
thin  airfoils  has  been  used.  The  blade  is  divided  into  NR  spanwise 
segments,  and  the  chordwise  aerodynamic  boundary  condition  is  satisfied 
at  the  midsection  of  each  of  these  segments  for  the  blade  in  each  of  NA 
azimuth  positions;  thus,  there  are  a  total  of  NR  A  (s  NR  *  NA)  points  in  the 
rotor  disk  where  the  boundary  condition  is  satisfied  and  the  airloads 
are  computed.  The  subscript  notation  for  these  points  over  the  disk  is 
given  in  Figure  3.  The  spanwise  loading  is  assumed  constant  over 
each  of  these  blade  segments,  and  it  is  from  the  ends  of  these  segments 
that  the  trailing  vorticity  of  the  wake  representation  streams  with 
strength  equal  to  the  change  in  strength  of  the  bound  vorticity  between 
adjacent  segments.  The  three-dimensional  effects  due  to  the  inter¬ 
action  of  a  spanwise  segment  with  its  neighbors  has  not  been  included. 
Thus,  each  blade  segment  is  considered  to  be  two-dimensional,  and 
the  three-dimensional  effects  of  the  wake  are  included  at  each  blade 
radial  station. 

The  chordwise  aerodynamic  boundary  condition  requires  that  the  sum 
of  the  velocities  normal  to  the  chord  be  zero  (i.  e.  ,  there  can  be  no 
flow  through  the  airfoil).  The  chordwise  distribution  of  these  normal 
velocities  is  composed  of  three  parts.  The  first  part  is  termed  the 
quasi-steady  part  and  is  due  to  orientation  of  the  blade  segment  relative 
to  the  total  local  velocity  vector,  its  camber,  and  its  pitching  and 
plunging  velocities  associated  with  the  blade  dynamic  response.  The 
second  part  is  termed  the  wake-induced  velocity  and  is  computed  from 
the  wake  vortex  representation.  The  third  part,  termed  the  airfoil- 
induced  velocity,  is  that  due  to  the  bound  vorticity  distribution  and  is 
required  to  balance  the  sum  of  the  first  two  parts.  That  is,  it  is  the 
bound  vorticity  distribution  which  is  considered  to  be  the  unknown  and 
must  be  such  as  to  create  an  airfoil-induced  velocity  distribution  which 
just  cancels  the  sum  of  the  quasi-steady  and  wake-induced  velocities. 

This  sum  of  the  quasi-steadv  and  wake-induced  velocities  is,  thus, 
termed  the  impressed  velocity  distribution  because,  when  it  is  impressed 
on  the  airfoil,  a  bound  vorticity  distribution  is  called  for  to  create 
canceling  airfoil-induced  velocities. 


The  coordinate  system  used  for  description  of  the  blade  representation 
is  shown  in  Figure  4  where  the  airfoil  chord  is  centered  on  the  %  -axis 
with  the  leading  and  trailing  edges  at  -b  and  + 1> ,  respectively.  For 
convenience  in  the  formulation,  a  transformation  is  made  to  a  secondary 
coordinate,  q  ,  by  the  following  transformation; 

pc  -  -  b  cos  &  ( 1 ) 
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Let  the  chordwise  distributicm  of  bound  vorticity  (i.  e.  ,  the  chordwise 
singularity  distribution  required  to  satisfy  the  boundary  condition)  on  a 
given  blade  radial  station,  for  each  instant  of  time,  be  expressed  in  a 
Glauert  series,  thusly: 

y(6,  -t)  -  Z  | An  (t)  cot  +  y  Antf)sin  n  0  J.  (2) 


The  total  circulation  on  the  airfoil  is  given  by 


Thus, 


/rb 

Y(%)cL% . 

-o 


r=  2rrb(A0+±At).  (3) 

The  chordwise  distribution  of  the  pressure  difference  is  given  by  the 
following  form  of  the  linearized  Bernoulli  equation  for  unsteady  flow, 

Zl p(*>t)  =/>[</,  Y(z,t)-b~  (4) 

and  can  be  expressed  in  terms  of  the  Glauert  coefficients  by  substituting 
(2)  into  (4),  thusly: 

Ap(G.t)  -  ZfiVt^AQcot^  +  £  An  sin  nO^ 

■f  2/3 +  z  At)0  +■  (Ao  +  j/\2  )  sin  9  (5) 


The  airfoil  unsteady  lift  and  pitching  moment  about  the  midchord  per 
unit  span  can  also  be  expressed  in  terms  of  Glauert  coefficients  by 
substituting  (5)  into  the  following  expressions: 


thus, 


Z*  J*b Ap(x)dx 

•r  +b 

77]  =  y  xApix)cLx,  i 


Z  =  Zrrbp  V,  \^A0  +  —A,  +  jpr  —  (3A0  +  Ai  +■  j  Az) j 
/>?  =  7Tbz/>  /,  \Ao+{Az-^Jd  (7) 
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It  should  be  noted  the  lift  depends  only  on  the  first  three  coefficients 
of  the  Glauert  series,  and  pitching  moment  depends  only  on  the  first 
four. 


Due  to  real  fluid  effects,  the  bound  circulation,  lift,  and  pitching 
moment  do  not  experimentally  attain  the  magnitudes  given  by  equations 
(3),  (6),  and  (7).  An  empirical  correction  factor,  Ci^jZrr  ,  is  there¬ 
fore  applied  to  the  bound  circulation  and  the  circulatory  part  of  the 
lift  and  moment;  the  Ci ^  is  the  measured  lift-curve  slope  for  the  air¬ 
foil  section.  The  bound  circulation,  lift,  and  pitching  moment  are  thus 


r  (8) 
£  =  bp  Vt  (A0  f~A ,)  +TrbpJ~  (3A0  +  Af+  j  A z  )  (9) 

^7  =  (Ao-h-A,}  +  -~7rb pVt  (-At  +A2) 

(10! 


The  computation  of  the  lift,  moment,  and  Glauert  coefficients  is  in  the 
time  domain  and  for  each  radial  station  the  Glauert  coefficients  are 
known  only  at  NA  equally  spaced  discrete  azimuth  positions.  Numerically 
the  time  derivatives  indicated  in  equations  (9)  and  (10)  were  evaluated 
at  these  discrete  points  by  a  procedure  which  "effectively"  fits  the 
points  with  a  Fourier  series  and  then  takes  the  slope  of  this  curve. 

This  procedure  is  presented  in  Appendix  III. 


Because  the  bound  vortex  strength  for  each  radial  segment  is  periodic 
ana  the  wake  vorticity  depends  on  the  radial  and  time  rates  of  change 
of  the  bound  vorticity,  each  of  the  strengths  of  the  wake  vortex  elements 
can  be  expressed  in  terms  of  the  NRA  bound  vortex  strengths  as  indicated 
in  Figure  5.  If  the  relative  orientations  of  the  blade  segments  and  the 
wake  elements  are  known  (or  prescribed),  application  of  the  Biot-Savart 
law  enables  the  chordwise  distribution  of  induced  velocity  on  a  blade 
segment  at  a  point,  K  ,  in  the  disk  to  be  expressed  in  terms  of  the  wake 
vortex  strengths  and,  thus,  in  terms  of  the  A/RA  strengths  of  the  bound 
vorticity  as 

*-*(*)  =Z  c«i (%)  o  (11) 

J-  f 

where  Ckj(%)  is  the  sum  of  the  Biot-Savart  coefficients  (expressing 
the  induced  velocity  distribution  over  the  chord  due  to  a  single  straight- 
line  element)  for  all  wake  elements  whose  strengths  depend  on  Q  . 

The  derivation  of  the  expression  for  the  chordwise-induced  velocity 
distribution  due  to  a  single  arbitrarily  oriented  vortex  filament  is 
presented  in  Appendix  II. 
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For  each  blade  segment  and  instant  of  time  (i.  e. ,  at  each  point,  K  ,  in 
the  rotor  disk),  the  expression  of  the  chordwise  aerodynamic  boundary 
condition  is  given  by  the  following  integral  equation: 


where 


---  f  — - 

17T  .  Xj  - 


-  (~h  +  V<  cc?  *  6c  z)K 


^(*)=ZcK.(x)r;  J 

J  “  ? 

are,  respectively,  the  quasi-steadv  and  wake-induced  parts  of  the 
chordwise  normal  velocity  distribution.  If  the  wake-induced  part  of  the 
normal  velocity  distribution  (expressed  in  the  G  -coordinate)  is  expanded 
in  a  Fourier  cosine  series 

HRA  eo 

-  (us 

j-t  »~1 

then  the  impressed  velocity  distribution  which  is  defined  as  the  sum 
of  the  quasi-steady  and  wake-induced  velocities  can  be  written  as 


where 


iyK(Q)  +  uxK(e)  =  +  Z  v'n  cos  nG 


^0K  -  (~h  +  Vi  ag)x  +  Z  0 


K  *  (-**,)«  +  E  •?<,.  rJ  <16> 

j’l  4 

MR* 

=  ZL  ^nK.  0  ’  n  ~ 

j  *  /■  J 

and  Snj..  is  the  sum  of  the  nth  harmonic  cosine  components  (of  the 
Biot-Savart  coefficients)  for  each  wake  element  whose  strength  depends 
on  F,-  .  The  derivation  of  the  expressions  for  the  first  four  cosine 
components  of  the  chordwise  normal  induced  velocity  distribution  (i.  e.  , 
<Sq  ,  Sf  ,  S2  and  Sj  )  due  to  a  single  wake  vortex  filament  is  presented 
in  Appendix  II. 

The  result  of  evaluating  the  integral  of  (12),  after  changing  to  the  &  - 
coordinate  and  substituting  (2)  for  y>((9)  is  also  a  cosine  series, 


t-L  f*b  Z*liM£  -it**,  ccne. 

*-«f  n?i 


II 


(18) 


By  substituting  (15)  and  (17)  into  (12), 

oo  oo 

V'ok  +H  *r.K  cosnd  =  A0k  A„k  cos  nQ, 

n‘*  ml 

and  equating  the  coefficients  of  like  terms,  the  Glauert  coefficients  can 
be  expressed  in  terms  of  the  coefficients  of  the  series  representing  the 
impressed  velocity  distribution,  thusly: 

//AG  ' 

Aoy  -  t^OK  =  Yf  +  2—-  ^°k-  O 

fJKA 

A<k  ~  {.+  b^s)x  "Z_,  Sh q  rj 

Jzf 

=  "  ~  ^rtK;  0 

J-t  J 


(19) 
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As  a  result  of  this  identification  of  Glauert  coefficients  with  the  cosine 
coefficients  of  the  impressed  velocity  distribution,  the  lift  and  moment 
are  seen  to  depend  on  only  the  first  four  coefficients  of  the  cosine 
expansion  of  this  impressed  velocity  distribution. 


From  (8)  and  (19),  the  equation  for  the  bound  vortex  strength  at  a 
point,  K  ,  of  the  disk  is 

,  HAA  . 

r*  -  (20) 

*  J  si  •> 

or  by  letting  .  j 


and 


°~Kj  -  CtocK  b*  (So*j  ~  2  S'kj  )  ? 
equation  (20)  becomes 


(21) 


r ;  I.  +• 


/=/  J 


(22) 


This  is  the  equation  for  the  bound  vortex  strength  of  each  blade  segment 
at  each  azimuth  position  (i.  e, ,  each  instant  of  time)  and,  thus,  the 
general  equation  of  a  simultaneous  set  of  NRA  equations  for  the  bound 
circulations.  The  constant  term,  IK  ,  of  each  of  these  equations 
depends  only  on  the  blade  response  and  flight  condition  (i.e.  ,  it  depends 
only  on  the  quasi-steady  part  of  the  impressed  velocity  distribution  and 
is  referred  to  as  the  quasi-steady  bound  vorticity);  the  second  term  is 
the  induced  part  of  /*-  .  This  set  of  equations  is  solved  iteratively 
using  the  Gauss-Sidel  method  (see  reference  4). 


For  the  present,  method,  the  only  blade  degrees  of  freedom  considered 
are  the  flapping  and  the  flapwise  bending;  the  pitching  and  torsion  motions 
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are  specified.  The  blade  normal  modes  are  used  in  the  description  of 
the  blade  dynamic  response  and,  thus,  the  total  plunging  displacement 
and  velocity  with  respect  to  the  rotor  shaft  plane  (plane  normal  to 
shaft  at  rotor  hub)  are 


HO 

4>K*(r)  r-.  JT  ts(r)qKA  (23) 

Szf  ’ 

Ar'C 

<P«»(r‘)  =  Z  fs(r)p*s  (24) 

5:/ 

where  MD  is  the  number  of  modes,  and  -f$(r )  and  qKAs  are,  respectively 
the  normalized  bending  mode  shapes  and  tip  deflections  in  each  mode,  5 
at  each  azimuth  station  KA  . 


The  h  term  in  the  I *  of  (21)  is  the  effective  plunging  velocity  and  includes, 
in  addition  to  the  plunging  velocity  of  (24),  the  component  of  the  rotor 
translational  velocity,  ,  normal  to  the  blade  at  each  spanwise 
segment,  thus: 

HO  HD 

h  =Vf  sin  OCS  +  (Vf  COS  Ocs)(cos  f/  (?)$KAs  +22*9  (f,)?KAs  -  (25) 

Szf  Szf 

When  (19)  is  substituted  into  (8),  the  equation  for  the  bound  vortex 
strength  is 

*K  =  2Lr°~’z  *  (26) 

The  linear  part  of  the  impressed  velocity  distribution  is  (vj,  +  v*,H  cos&). 
Therefore,  by  use  of  (1),  the  factor  (ya  -  ~  irf)  of  (26)  can  be  interpreted 
as  the  vaL  e  of  the  linear  part  of  the  impressed  velocity  distribution  at 
*  =  +  -§■  the  three-quarter  chord;  that  is,  the  bound  vorticity  depends  only 
on  the  linear  component  of  the  impressed  velocity  distribution  at  the 
three-quarter  chord.  This  fact  is  used  in  defining  the  airfoil  section 
stall  in  the  computation.  A  blade  section  is  defined  to  be  stalled  when 
the  linear  component  of  the  impressed  velocity  distribution  at  the  three- 
quarter  chord  exceeds  that  on  the  airfoil  at  its  measured  stall  angle  of 
attack.  Thus,  if  a ^  is  the  measured  stall  angle  for  the  section,  then 
the  value  of  the  bound  vorticity  at  stall,  is  by  (26),  rm*  =  Cj  bK  (/,  sinam) 
and  the  section  is  said  to  be  stalled  when  ' 

"  T  ~  J  V'lK 

- - — — — - - -  =  — - : -  >  1.0  .  (27) 

C^aKbt((y<KSln  m )  V,K  Stn  a m 


The  rmK  can  be  computed  once  the  value  of  <xm  is  specified.  During  the 
iterative  solution  of  the  simultaneous  set  of  r  -equations,  the  values  of 
the  rK  are  not  allowed  to  exceed  their  respective  f’m  .  Thus,  the 

K 
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magnitudes  of  the  r  which  are  shed  into  the  wake  are  limited.  This  is 
an  approximation  to  prevent  large,  unrealistic  values  of  vorticitv  from 
being  shed  into  the  wake  from  sections  which  are  above  their  stalling 
angle  of  attack* 

The  lift  and  pitching  moment  are  computed  by  (9)  and  (10)  when  R  <  1.0 
(i.  e. ,  when  che  airfoil  section  is  below  stall)  but,  when  R  >  1.0  and  the 
section  is  stalled,  the  lift  and  moment  are  approximated  by  the 
following: 


*/>v,K  r„K  +&JK  (28) 

where 

AU,k=  Z/>bK  VTk  sin 

and 

irix  =  ~kK£K  .  (29) 

That  is,  the  stalled  lift  is  computed  as  the  sum  of  the  stalled  (limited) 

circulatory  lift  and  a  cross -flow  drag  force,  Del  .  The  stalled  pitching 
moment  about  the  midchord  is  computed  as  the  moment  due  to  the  stalled 
lift,  ,  acting  at  the  quarter -chord.  The  velocity,  VT  ,  is  the  total 
local  resultant  velocity  acting  normal  to  the  spanwise  axis  (or  midchord 
axis)  of  the  local  blade  segment  and  is  given  by 

Vj  --  VfZ  +  VzZ  (30) 

The  velocities,  V,  and  V2  are  the  two  orthogonal  components  of  Vr  in  a 
plane  normal  to  the  local  spanwise  axis.  Velocity  K/  lies  on  the  inter¬ 
section  of  this  plane  with  a  plane  normal  to  the  shaft. 

Vt  =  (Vf  cos  a.$  )  sin  if)  +  Sir  (31 ) 

'Z£  =  ur  -  Vf  sin  a:s  -  {Vf  cosocs)(-~;)  cos  f)  -  j>  (32) 

where  <fi  is  the  blade  displacement  given  by  (23). 

The  angle,  <xe  (effective  angle  of  attack),  is  the  local  angle  between  the 
line  of  zero  lift  for  the  airfoil  and  the  direction  of  VT  ;  it  is  given  by 

cxe  =  +  &  (33) 

where  S  =  arctan  (V2/Vt )  ■ 
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COMPUTATIONAL  PROCEDURE 


The  computational  procedure  is  fundamentally  an  iterative  procedure 
which  repetitively  computes  the  aerodynamic  loads  and  solves  the  equa¬ 
tions  of  motion  of  the  blade  so  as  to  obtain  consistent  airloads  and  blade 
response.  The  resulting  solution  includes  all  the  aerodynamic  coupling 
between  the  degrees  of  freedom. 

The  computational  procedure  requires  that  the  wake  configuration  be 
specified;  that  is,  that  the  relative  orientation  of  the  blade  segments 
and  the  wake  vortex  elements  be  known.  To  establish  this  relative 
orientation,  the  first-order  displacements  of  the  blades  relative  to  the 
shaft  are  needed.  Thus,  in  addition  to  the  rotor  translational  velocity, 
shaft  angle,  and  wake  transport  velocities,  a  reasonable  approximation 
of  the  steady-state  displacement  of  the  blade  relative  to  the  shaft  (i.  e.  , 
coning  and/or  steady  bending),  and  the  first  harmonic  cosine  flapping 
angle  are  needed  to  establish  the  relative  orientation  of  the  blade 
segments  and  wake  elements.  If  the  steady  blade  displacement  and 
first  cosine  harmonic  flapping  (relative  to  the  rotor  shaft)  predicted 
by  the  method  are  significantly  different  from  the  estimates  used  to 
establish  the  relative  orientation  of  the  blades  and  wake,  the  computed 
values  should  be  used  to  re-establish  the  wake  orientation  and  then  the 
computation  repeated. 

First,  the  blade  and  wake  coordinates  are  computed,  and  then  an 
integration  (summation)  over  the  wake  is  performed  to  obtain  the  matrices 
of  induced  velocity  coefficients  SOK.  through  S3K.  as  defined  in  equation 
(14).  These  coefficients  are  for  the  wake-induced  velocities  normal  to 
the  tip -path-plane,  but  are  computed  at  the  steady-state  deflected 
positions  of  the  blades.  Because  the  blade  spanwise  slopes  are  small, 
these  induced  velocities  are  assumed  to  act  in  the  direction  of  V2  [  equation 
(32)].  From  these  induced  velocity  coefficients,  the  matrix  of  coeffi¬ 
cients,  ,  for  the  simultaneous  set  of  bound  vortex  equations  (  T  - 

equations)  is  computed. 

Next,  an  iterative  procedure  is  used  to  solve  the  equations  for  the  aero¬ 
dynamic  loads  and  blade  response.  (A  flow  diagram  indicating  the  main 
conceptual  elements  of  this  iterative  procedure  is  presented  in  Figure  6.) 
The  procedure  begins  with  a  first  approximation  to  the  blade  response 
and  then  evaluates  the  constants,  I  ,  [given  by  equation  (21)]  of  the  P  - 
equations.  The  P  -  equations  are  next  solved  iteratively  by  use  of  the 
Gauss-Sidel  method  and  then  the  Glauert  coefficients  evaluated  by 
equation  (19).  From  the  Glauert  coefficients,  the  lift  ana  pitching 
moments  are  computed,  and  then  a  "generalized  airload"  is  computed  in 
each  mode  for  each  azimuth  position  by  integrating  the  product  of  mode 
shape  and  the  spanwise  airload  over  the  span.  The  generalized  airload 
in  each  mode  is  then  harmonically  analyzed  and  used  as  the  forcing 
function  for  that  mode  at  each  harmonic.  The  blade  equations  of  motion 
for  each  mode  are  then  solved  for  the  response  at  each  harmonic.  With 
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this  new  approximation  to  the  blade  response,  the  Is  of  the  r  -  equa¬ 
tions  are  recomputed  and  the  entire  computational  process  is  repeated. 
This  process  is  repeated  until  the  blade  response  converges  to  within 
prescribed  limits. 

The  generalized  airload  is  not  the  forcing  function  as  conventionally- 
defined  because  it  contains,  in  addition  to  the  aerodynamic  forces  inde¬ 
pendent  of  the  response,  all  the  aerodynamic  loads  which  are  directly 
and  indirectly  a  result  of  the  blade  response.  These  additional  loads  are 
the  aerodynamic  spring,  damping  and  mass  forces,  and  the  aerodynamic 
coupling  force  from  the  wake.  Normally,  terms  representing  these  aero¬ 
dynamic  forces  that  are  proportional  to  the  system  reponse  are  put  on 
the  left  side  of  the  equation  of  motion  along  with  the  other  terms  which 
are  also  proportional  to  the  system  response.  The  forces  that  are 
independent  of  the  system  response  are  put  on  the  right  side  of  the  equa¬ 
tion.  In  the  computational  procedure  which  has  been  developed,  these 
response-dependent  components  of  the  airload  are  implicit  in  the  solution 
of  the  r  -  equations  and  the  computation  of  the  resultant  airloads.  There¬ 
fore,  by  using  the  generalized  airload  as  the  forcing  function,  the 
response-dependent  airloads  are  effectively  taken  from  the  left  side  of  the 
equation  and  included  with  the  conventional  forcing  function  on  the  right. 
This,  of  course,  does  not  alter  the  equation  of  motion  nor  its  solution. 

In  the  computational  procedure  described  above,  the  generalized  airload 
(right  side  of  equation  of  motions}  is  evaluated  for  the  existing  approxi¬ 
mation  to  the  response  (from  previous  iteration).  By  using  this  as  the 
forcing  funct  on,  the  equation  of  motion  is  then  solved  for  the  next 
approximation  to  the  response.  Unfortunately,  when  this  process  is 
repeated  it  will  not  always  converge. 

An  investigation  was  made  to  determine  and  understand  the  controlling 
factors  of  this  iterative  procedure.  It  was  found  that  convergence  of  this 
procedure  was  improved  when  the  dependence  of  the  right  side  of  the 
equation  on  the  variables  (i.  e.  ,  on  the  response)  was  reduced.  This  is 
evident  when  the  right  side  is,  in  fact,  independent  of  the  variables,  then 
the  result  of  the  first  iteration  is  the  solution.  Thus,  in  the  computational 
procedure  which  has  been  developed,  the  dependence  of  the  right  side 
(i.  e.  ,  the  generalized  airloads)  on  the  response  is  reduced  by  subtracting 
an  approximation  for  the  quasi -steady  damping  from  both  sides  of  the 
equation.  That  is,  a  term  which  is  believed  to  be  the  largest  response- 
dependent  part  of  the  airload  has  effectively  been  brought  back  to  the  left 
side  of  the  equation  from  the  right  side.  This  quasi-steady  damping 
subtracted  from  the  right  side  is  computed  as  a  function  of  the  response 
of  the  previous  iteration,  while  that  subtracted  from  the  left  is  computed 
as  a  function  of  the  response  of  the  present  iteration.  When  the  solution 
converges,  these  two  are  equal  and  the  equation  is  unaltered.  This 
procedure  has  worked  very  well,  thus  far,  and  the  solution  obtained  is 
independent  of  this  approximate  quasi-steady  damping  correction  term, 
but  the  rate  of  convergence  does  depend  on  it.  It  was  found  that  both 
halving  and  doubling  this  term  increased  the  number  of  iterations  required 
for  convergence  without  changing  the  solution. 
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The  equations  of  motion  for  the  blade  which  were  actually  used  in  the 
.iterative  procedure  are  obtained  as  indicated  in  the  following.  The 
equation  of  motion  for  a  generalized  coordinate  <gs  (tip  deflection  in  a 
normal  mode)  expressing  the  balance  of  all  the  generalized  forces  in 
that  coordinate  is 

Ms  fs  +  ‘S,  is  (34) 

where  the  first  term  on  the  left  is  the  inertia  force,  the  second  is  the 
structural  damping  force  (proportional  to  the  spring  force  but  in  phase 
with  the  velocity,  os  is  the  structural  damping  coefficient),  the  third 
is  the  spring  force,  and  the  term  on  the  right  is  the  generalized  ai-load. 
The  quasi-steady  damping  term  to  be  subtracted  from  both  sides  of 
equation  (34)  to  improve  convergence  is  approximated  as 

Ds  ~  ~  J*  b(r)rf*(r)  dr  (35) 

The  generalized  coordinate  and  the  generalized  airload  are  expanded  in 
the  following  Fourier  series: 

oo 

-  a>0^  t  (q>n^  COS  n  Hi  +  sir?  nSXi  )  (36) 


Gc  =  An  t 


y*  cos  n  Cli  -f-  B n  Sinn  Sit) 


By  substituting  these  expressions  for  <fs  and  Gs  in  equations  (34)  and 
(35),  subtracting  equation  (35)  from  both  sides  of  equation  (34),  and  then 
equating  coefficients  of  like  terms  in  the  result,  the  following  are  obtained 
as  the  equations  of  motion  used  in  the  iterative  procedure: 


.  -  Z  d)  .  ( t-f ) 

/%  oos  0,0  -  A0 

.  .  /  Z  2  _  2  \  (C) 


K(u>s~n  a  )a'n!  9s  +  Czr'SX)bn,)~Aol 


Js.*  ~  &n  f)+(Csnn.)o!yV; 

n-  i,z,  3.  •  -• 

where  the  superscript  (£)  indicates  the  iteration  number  and 


r*  , 

Cs=Cj  /SJXf  rb(r)fs  (r)ctr 
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is  the  coefficient  of  the  quasi-steady  damping  term.  Thus  for  each 
generalized  coordinate  (degree  of  freedom)  at  each  harmonic  of  the 
rotational  speed  there  is  a  pair  of  simultaneous  equations  for  the  sine 
and  cosine  components  of  the  response. 

After  the  iterative  procedure  has  converged,  the  radial  distribution  of 
the  structural  bending  moments  at  each  harmonic  (cosine  and  sine)  are 
computed  from  the  blade  response  in  each  normal  mode  and  the  normal 
mode  moment  distributions  as  follows: 

NO 

M0(r)  =  C'o  ff)s(r) 
s*t  S 
_  NO 

Jr)  =  A.  ma(r) 

s=f  s 
UD 

=  Z  bns™s (r) 

s-f  5 

where  the  ms(r)  are  the  bending  moment  distributions  in  each  normal 
mode,  &  .  It  should  be  noted  that  in  this  procedure  the  steady  blade 
displacements  are  represented  by  the  sum  of  the  steady  deflections  in 
the  normal  modes.  While  the  number,  A ID ,  of  modes  used  is  adequate 
to  represent  the  significant  dynamic  response  of  the  blade  (because  the  , 
frequency  response  in  each  mode  is  such  that  significant  dynamic 
response  occurs  only  near  the  natural  frequency  of  the  mode),  this  l 

number  of  modes  is  not  necessarily  adequate  to  represent  the  static 
deflection  shape  of  the  blade  for  purposes  of  computing  steady  bending 
moments.  This  procedure  will,  however,  yield  the  static  blade  deflec¬ 
tion  shape  with  sufficient  accuracy  for  the  aerodynamic  aspects  of  the 
problem.  This  is  because  the  blade  section  effective  angle  of  attack  is 
a  function  of  the  blade  spanwise  slope  (first  derivative  with  respect  to 
radius)  whereas  the  bending  moments  depend  on  the  local  blade  curvature 
(second  derivative  with  respect  to  radius).  For  a  more  accurate  pre¬ 
diction  of  the  steady  bending  moments,  corresponding  to  the  steady-state 
aerodynamic  loads  (plus  the  centrifugal  loading  due  to  built-in  precone, 
if  used),  one  of  the  standard  procedures  for  computing  static  deflections 
of  beams  should  be  used. 
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COMPUTER  PROGRAM 


The  computer  program  is  written  entirely  in  Fortran  IV  for  the  IBM 
7044  EDP,  and  is  physically  two  separate  programs  (Part  1  and  Part  2) 
which  are  run  sequentially.  Each  part  required  essentially  the  entire 
available  core  storage  capacity.  The  program  required  two  tape  drives 
in  addition  to  the  three  for  the  monitor  system  input,  output,  and  library 
files.  Each  of  these  two  parts  is  subdivided  into  subroutines  for  ease 
and  convenience  of  programming,  check  out,  and  the  incorporation  of 
modifications.  The  running  time  for  the  complete  program  depends  on 
the  individual  problem  being  processed  (i.  e.  ,  on  the  number  of  blades, 
blade  segments,  azimuth  positions,  the  wake  configuration  used,  the 
number  of  revolutions  of  wake,  the  number  of  degrees  of  freedom  used, 
and  the  number  of  iterations  required  to  obtain  the  solution).  For  the 
computations  which  have  thus  far  been  made;  Part  1  of  the  program  has 
required  from  0.  10  to  0.  20  hour  and  Part  2  has  required  from  0.  15  to 
0.  25  hour  of  machine  time. 

Part  1  of  the  program  generates  the  blade  segment  and  wake  element 
coordinates,  and  then  it  evaluates  the  matrices  of  induced  velocity 
coefficients  by  integrating  the  Biot-Savart  expression  over  the  wake  — 
these  are  stored  on  tape.  Part  2  of  the  program  evaluates  the  coeffi¬ 
cients  for  the  r  -equations  and  the  equations  of  motion;  then  it  solves 
the  system  of  equations  by  the  iterative  procedure.  After  the  solution 
has  converged,  the  time  histories  of  all  the  results  are  harmonically 
analyzed  and  both  the  time  histories  and  their  Fourier  coefficients  are 
printed  as  outputs  together  with  all  the  inputs  for  both  parts  of  the 
program.  Flow  diagrams  indicating  the  logic  for  both  parts  of  the 
program  are  presented  in  Appendix  IV, 
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COMPUTED  RESULTS  AND  COMPARISONS 

-  WiTHMEA'SuKEIJ'R'ESULTS - 


The  computational  procedure  was  used  to  analyze  two  flight  conditions  of 
each  of  two  rotor  systems  for  which  measured  airloads  and  blade  bend¬ 
ing  moments  are  available.  These  four  computations  are  for  the  HU-1A 
rotor  at  advance  ratios  of  0.  08  and  0.  26,  and  for  the  H-34  rotor  at 
advance  ratios  of  0.  18  and  0.  29.  Comparisons  of  the  measured  and 
computed  lift  loads  and  blade-flapwise  bending  moments  are  presented; 
however,  only  the  computed  pitching  moments  are  presented  because 
measurements  of  them  are  not  available. 

The  dynamic  response  of  a  rotor  blade  can  be  expressed  in  terms  of  its 
response  in  each  of  its  normal  modes.  The  response  in  each  normal 
mode  at  each  harmonic  of  the  rotational  speed  is  proportional  to  both  the 
amplification  factor  and  the  generalized  force  for  the  specific  mode  and 
harmonic.  The  amplification  factor  is  determined  by  the  ratio  of  the 
harmonic  forcing  frequency  to  the  natural  frequency  (of  the  mode  under 
consideration),  and  the  generalized  force  is  determined  by  the  mode 
shape  and  radial  distribution  of  the  airload  at  the  specific  harmonic. 
Therefore,  to  predict  accurately  the  bending  moments  at  each  harmonic 
(and,  thus,  the  total  bending  moments),  it  is  necessary  to  know  the 
natural  frequencies  and  to  predict  the  radial  distribution  of  the  airload 
at  each  harmonic  accurately;  it  is  not  sufficient  to  simply  predict  time 
histories  of  the  airload  which  look  similar  (that  is,  agree  approximately 
in  the  lower  harmonics)  to  the  measured  airloads  at  a  few  radial  stations. 
Because  of  the  importance  of  the  harmonic  components  of  the  airload, 
their  radial  distributions  are  presented  here,  in  the  comparisons  of  the 
measured  and  computed  results,  in  addition  to  the  time  histories  of  the 
total  airload.  Furthermore,  presentation  of  the  harmonics  is  a  much 
more  detailed  view  of  the  airload  distribution  (over  the  rotor  disk) 
because  it  displays  its  constituents  and  can  thus  provide  a  much  more 
critical  comparison.  Both  forms  of  presentation  are  also  used  for  the 
bending  moments. 

In  addition  to  the  computations  made  for  comparison  with  the  measured 
results,  a  computational  investigation  was  made  to  determine  the 
relative  sensitivity  to  some  of  the  input  wake  configuration  parameters. 

It  v/as  found  that  the  sensitivities  depend  on  the  particular  application. 

The  only  general  conclusion  which  can  be  drawn  is  that  the  wake  configur¬ 
ation  is  important  and  a  better  means  for  establishing  it  more  accurately 
is  needed.  Results  of  this  investigation  are  presented  later. 

All  the  input  information  required  for  the  computations,  as  well  as  the 
measured  results  used  for  the  comparisons,  was  obtained  from  refer¬ 
ences  5  and  6.  The  bending  normal  mode  shapes,  moment  distributions, 
and  natural  frequencies  were  computed  from  the  blade  mass  and  elastic 
distributions  contained  in  these  references.  All  the  computations  for 
both  rotor  systems  were  made  for  24  blade  azimuth  positions.  For  the 
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computations,  the  HU-1A  blade  was  divided  into  8  spanwise  segments 
and  the  H-34  into  9  spanwise  segments  as  shown  in  Figure  7.  The 
midpoints  of  these  segments  include  all  the  radial  stations  for  which 
measured  airloads  are  available.  In  all  the  computations  presented  for 
comparison  with  measured  results,  the  first  three  azimuthal  segments 
of  wake  behind  each  blade  were  of  the  grid  representation  and  the  remain¬ 
der  of  the  wake  was  only  a  tip  vortex  positioned  at  90  percent  (represent¬ 
ative  contraction  observed  in  reference  2)  of  the  blade  radius.  The 
number  of  revolutions  of  wake  varied  for  each  case  depending  on  the 
advance  ratio.  Unless  otherwise  noted,  the  computed  momentum  induced 
velocity  is  used  for  the  induced  part,  ur(r,  f)  ,  of  the  Z  -component  of 
wake  transport  velocity. 

The  steady  components  (mean  value)  have  been  subtracted  from  the  time- 
history  presentations  in  the  comparisons  of  computed  and  measured 
results,  and  the  radial  distribution  of  these  components  is  presented 
with  the  other  harmonic  components. 
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HU-1A  at  -  0.  26 


To  establish  the  wake  configuration  for  this  case,  the  following  are  the 
values  of  the  parameters  -used: 

Vf  =  188.0  feet  per  second 

CCj  -  5. 8  degrees 

ur(r , tp)  -  -6.0  feet  per  second 

NW  =  2  (Number  of  revolutions  of  wake) 

ADV  =  0.  7  (wake  advance) 


The  measured  and  computed  time  histories  and  radial  distributions  of  the 
harmonics  of  the  lift  loadings  and  the  blade-flapwise  bending  moments 
are  presented  in  Figure  8.  From  the  time  histories,  it  appears  that 
there  is  reasonable  correlation  between  the  measured  and  computed 
results.  The  correlation  does  not,  however,  appear  as  good  when  the 
radial  distribution  of  the  harmonics  is  compared.  In  general,  the 
correlation  for  the  airloads  is  best  for  the  lower  harmonics.  The  radial 
distributions  of  the  higher  harmonics  of  the  airload  are  not  as  well 
predicted  as  their  relative  amplitudes.  This  is  true  at  most  harmonics 
of  the  bending  moments.  It  is  noted  that  above  the  5th  harmonic  the 
measured  and  computed  bending  moments  are  insignificant. 

The  reversal,  at  the  inboard  end  of  the  blade,  in  the  radial  distribution 
of  the  measured  steady-bending  moment  is  due  to  the  centrifugal  force 
moments  (arising  from  the  built-in  preconing  of  the  blades)  exceeding  the 
lift  moments  at  the  inboard  end  of  the  blade.  Because  the  steady 
component  of  the  lift  loading  is  overpredicted,  the  lift  moment  exceeds 
these  centrifugal  moments  at  the  inboard  end  of  the  blade  and  the  reversal 
does  not  occur  in  the  computed  steady-bending  moment. 

In  the  time  history  of  the  airloads,  it  was  observed  that  the  computed 
airloads  are  slightly  advanced  with  respect  to  the  measured  airloads  on 
the  advancing  side  of  the  disk  and  retarded  on  the  retreating  side.  This 
effect  would'  result  if  the  wake  configuration,  used  in  the  computation, 
was  skewed  slightly  more  than  the  actual  wake.  Thus,  a  blade  would 
see  the  wake  of  the  preceding  blade  sooner  on  the  advancing  side  and 
later  on  the  retreating  side.  This  small  difference  could  arise  because 
of  an  error  in  the  advance  ratio  or  because  of  induced  velocities  in  the 
X  -direction  (which  have  been  neglected),  increasing  the  transport 
velocity  in  that  direction. 
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HU  - 1 A  at  fic  -  0.  08 


The  wake  configuration  for  this  case  was  established  by  the  following 
values  of  the  parameters: 

Vf  -  55.  1  feet  per  second 

<X7  -  2.  5  degrees 

cv(r, p)  -  -21.0  feet  per  second 

NY/  -  6  (revolutions  of  wake) 

AW  -  0.7  (wake  advance) 


The  time  histories  and  radial  distribution  of  the  harmonics  of  the  air¬ 
loading  and  blade -flapwise  bending  moments  are  presented  in  Figure  9. 
The  correlation  here  is  not  as  good  as  that  obtained  at  /^  =  0.  26.  The 
mean  value  of  the  airloading  has  been  underpredicted  and  likewise  the 
higher  harmonics  of  the  airloading  have  been  underpredicted.  The 
discrepancies  are  attributed  primarily  to  the  wake  configuration's  being 
a  poor  representation  of  the  real  wake  for  this  case.  The  abrupt 
changes  in  the  airload  on  both  the  advancing  and  retreating  sides  of  the 
disk  (which  are  due  to  the  induced  velocity  variation)  are  not  well  pre¬ 
dicted  for  this  case. 


The  various  contributions  to  the  computed  lift  loading  at  r/R  =  0.  95  are 
presented  in  Figure  10.  Here  in  Figure  10a  are  presented,  for  this 
radial  station,  the  quasi-steady  bound  vorticity,  J  ,  and  the  total 
bound  vorticity,  F  ,  together  with  the  lift  loading  X  .  The  induced 
velocity  (which  is  directly  proportional  to  the  difference  between  F  and 
the  quasi-steady  bound  vorticity)  is  presented  in  Figure  10b.  It  is 
observed  here  that  the  induced  velocity  is  most  negative  over  the  rear 
half  of  the  disk  and  is  nearly  zero  or  positive  over  the  forward  half  of 
disk.  This  is  characteristic  of  the  induced  velocity  distributions  of 
translating  rotors  because  (as  observed  in  Figure  21)  the  blade  stations 
are  generally  on  the  downwash  side  of  the  preceding  blade  tip  vortices 
as  the  blade  goes  from  f  -  270  degrees  to  V  -  90  degrees  (rear  half  of 
disk)  and  on  the  upwash  side  of  them  as  the  blade  goes  from  V  =90 
degrees  to  V  -  270  degress  (forward  half  of  disk).  It  is  also  observed 
that,  while  the  induced  velocity  incurs  a  large  total  change  as  the  blade 
passes  from  the  rear  half  of  the  disk  to  the  forward  half  of  the  disk,  it 
does  not  change  as  abruptly  as  would  be  expected  (nor  as  is  observed  in 
Figure  13b,  for  the  H-34  at  /x.  =  0.  18).  The  real  wake  (tip  vortex)  is 
probably  much  closer  to  the  rotor  disk  on  the  advancing  and  retreating 
sides  of  the  disk  because  of  the  Z  -  component  of  the  distortion  of  the  tip 
vortex  due  to  its  rolling  up  on  itself  at  the  sides  of  the  wake  as  described 
in  the  section  on  the  wake  representation.  Also,  corresponding  distor¬ 
tions  in  the  X  and  Y  -  directions  could  make  the  approach  of  the  blade 
to  the  preceding  blade's  tip  vortex  more  abrupt  (i.  e. ,  the  angle  between 
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this  vortex  and  the  blade  would  be  small)  as  viewed  in  X~Y  plane  rather 
than  very  gradual  as  was  the  case  for  this  computation  (where  the  angle 
is  nearly  90  degrees).  The  reversal  in  the  distribution  of  the  steady¬ 
bending  moment.  Figure  9d,  and  negative  value  near  the  blade  root  are 
overpredicted  for  this  case  because  the  steady  component  of  lift  was 
underpredicted  and  thus  the  negative  moments  due  to  the  centrifugal 
loading  arising  from  the  preconing  are  predominant.  The  overall  correl¬ 
ation  of  the  harmonics  of  the  bending  moments  appears  to  be  fairly  good 
compared  to  the  correlation  of  the  airloads.  This  is  because  the  higher 
harmonics  of  the  bending  moments  are  negligible,  and  it  is  the  higher 
harmonics  of  the  airload  which  are  not  well  predicted.  The  correlation 
of  the  lower  harmonics  of  the  airload  is  reasonable  and  so  also  are  the 
significant  lower  harmonics  of  the  bending  moments. 


H-34  at  /6  =  0.  29 

The  following  parameter  values  were  used  to  establish  the  wake  configur¬ 
ation  for  this  case: 

Vf  =  189.  2  feet  per  second 

<zT  -  6. 0  degrees 

ur(W)  -  -6.4  feet  per  second 

a/W  =  2  (revolutions  of  wake) 

The  comparisons  of  the  measured  and  computed  time  histories  and 
harmonic  distributions  of  the  airloads  and  bending  moments  are  pre¬ 
sented  in  Figure  11.  Because  of  errors  in  the  measured  steady  compon¬ 
ent  (zeroth  harmonic)  of  the  bending  moments,  they  are  not  presented  in 
Figure  lid.  The  time  histories  of  the  airload  appear  to  correlate 
reasonably  well.  However,  the  harmonic  distributions  indicate  that  the 
correlation  is  primarily  at  the  lower  harmonics,  while  the  higher 
harmonics  are  generally  overpredicted,  it  is  interesting  to  note  that 
the  relative  correlation  of  the  first  harmonic  cosine  component  of  the 
lift  is  very  similar  to  that  obtained  for  the  HU-1A  at  ju,  =  0.  26,  and  that 
the  first  and  second  harmonic  components  of  the  lift  for  these  two 
cases  are  also  similar. 

It  is  evident  from  the  time  histories  of  the  bending  moments,  that  only 
the  low  harmonics  correlate,  and  this  is  confirmed  in  the  harmonic 
distributions.  The  higher  harmonics  are  generally  overpredicted.  This 
is  especially  true  for  the  eighth  harmonic  and  is  due  to  the  third  bending 
mode  natural  frequency  being  near  the  eighth  harmonic  of  the  rotational 
speed. 
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H-34  at  =0.18 


For  this  case,  the  wake  configuration  was  established  by  the  following 
values  of  the  wake  parameters: 


aT 


ur(r,  (p) 
NW 
ADV 


118.2  feet  per  second 
3. 3  degrees 
-10.  6  feet  per  second 

3  (revolutions  of  wake) 
0.7  (wake  advance) 


The  time  histories  and  the  radial  distribution  of  harmonics  of  the  lift 
loading  and  bending  moments  are  presented  in  Figure  12.  The  correl¬ 
ation  in  both  the  time  histories  and  the  harmonics  of  the  airloads  for  this 
case  is,  in  general,  better  than  for  the  other  cases.  The  abrupt  rise 
and  fall  of  the  airloads  on  the  advancing  and  retreating  side  of  the  disk, 
Figure  12a,  have  apparently  been  slightly  overpredicted  here,  whereas 
for  the  low  advance  ratio  case  for  the  HU-1A  they  were  underpredicted. 

It  is  also  noted  that  the  computed  airloads  are  slightly  advanced  with 
respect  to  the  measured  airloads  on  the  advancing  side  of  the  disk  and 
slightly  retarded  on  the  retreating  side  as  was  noted  for  the  HU-1A  at 
/w.  =  0.  26.  This  seems  to  be  generally  true  for  most  of  the  computations. 
Here  again  the  harmonic  distributions  portray  the  discrepancies  between 
the  measured  and  c-  mputed  airloads  while  the  time  histories  appear 
to  correlate  reasonably  well.  The  computea  radial  distributions  of  the 
steady  and  first  harmonic  components  compare  quite  well  with  those 
measured.  The  computed  radial  distributions  of  the  "econd.  and  higher 
harmonics  contain  relatively  large  radial  variations.  Similar  varia¬ 
tions  are  evident  in  the  second  and  third  measured  components  but  with 
a  small  amplitude. 


Because  of  errors  in  the  measured  »teady  component  (zeroth  harmonic) 
of  the  bending  moments,  they  are  not  presented  in  Figure  12d.  The 
oscillatory  moments  correlate  reasonably  well,  except  for  the  higher 
harmonic  response,  due  primarily  to  the  third  bending  mode. 

An  interesting  insight  into  the  genesis  of  the  airloading  can  be  obtained 
from  Figure  13,  where  for  r/R  =  0.  95,  the  quasi-steady  component,  J  , 
of  the  bound  vorticity,  the  total  bound  vorticity,  r  ,  the  airloading,  £  , 
and  the  wake-induced  velocity  time  histories  are  presented.  It  should 
be  remembered  that  the  difference  between  the  I  and  P  curves  is  the 
induced  part  of  the  bound  vorticity  which  is  proportional  to  the  wake- 
induced  velocity.  First,  it  is  noted  for  this  case  that  the  decreasing 
airload  as  the  blade  goes  from  if/  -  0  to  ^  =  60  degrees  is  not  due  to  the 
blade  encountering  abruptly  the  downwash  of  a  trailing  vortex,  but  is 
due  to  the  rapidly  decreasing  quasi-steady  component  of  P  (i,  e,  ,  due 
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to  the  blade  velocity,  motions,  and  pitch  inputs).  The  sharp  rise  in 
the  airload  as  the  blade  goes  from  =  300  degrees  to  ^  =  330  degrees 
is  similarly  due  to  I  .  The  sharp  decrease  from  =  75  degrees  to 
f  =  120  degrees  and  the  rise  from  V  =  240  to  V  =  290  degrees  are 
partially  due  to  the  wake-induced  velocity  variation,  but  predominately 
the  result  of  the  I  variation.  The  sharp  rise  on  the  advancing  side, 
and  drop  on  the  retreating  side,  is  due  only  to  the  wake-induced  velocity 
variation.  Thus,  it  is  seen  that  the  character  of  the  airloads  is  strongly 
influenced  by  the  quasi-steady  aspect  of  the  problem.  The  wake-induced 
velocity  variations  are  characteristically  "hat  shaped",  as  in  Figures 
10b  and  13b. 

In  the  past,  there  has  been  speculation  that  the  lift  variation  on  the 
advancing  and  retreating  sides  of  the  disk  (which  resembles  the 
"characteristic"  induced  velocity  variation  of  a  single  vortex)  was  a 
reflection  of  the  induced  velocity  field  of  the  tip  vortex  of  the  preceding 
blade.  From  the  above  discussion,  it  is  evident  that  this  is  not  necessar¬ 
ily  true.  For  this  case,  the  lift  variation  (which  has  this  characteristic) 
is  the  result  of  both  the  wake-induced  and  quasi-steady  velocities  and 
not  primarily  the  reflection  of  only  the  induced  velocity  experienced  by 
the  blade  section. 


Pitching  Moments 


The  measured  pitching  moments  were  not  available;  therefore,  only 
the  computed  pitching  moment  time  histories  are  presented  in  Figure 
14.  It  should  be  noted  that  these  pitching  moments  are  about  the  airfoil 
midchord.  The  computed  pitching-moment  time  histories  are  seen  to 
be  very  similar  to  the  corresponding  lift  time  histories;  this  is  because 
the  predominant  part  of  the  pitching  moment  about  the  midchord  is  due 
to  the  lift  acting  at  the  quarter-chord.  This  can  be  seen  by  evaluating 
the  pitching  moments  about  the  quarter -chord  X--b/2\  thus, 


Mr 


b  r*  b 

*~zs~Jb  (*  +  %)Ap(x)dx 

XAp(x)  dx  -  ~J  A  p  (x)  dx  • 


(41) 


From  equations  (6)  and  (7),  these  two  terms  are  identified  as  the  pitch¬ 
ing  moment  about  the  midchord,  M  ,  and  the  moment  of  the  lift,  £  , 
acting  at  the  quarter -chord,  thus  the  moment  about  the  midchord  is 

M-.  +  \X.  (42) 


The  moment  about  the  quarter-chord  can  be  expressed  in  terms  of  the 
Glauert  coefficients  by  substituting  equation  (5)  into  (41)  to  obtain 
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^bZfiVM^A2)-Tb3pJi  [«„  H3) 

In  Figure  15,  for  r/R  -  0.75  and  0.  95  on  the  H-34,  the  computed  pitch¬ 
ing  moments  about  the  midchord  are  compared  with  the  moment  due  to 
the  lift  acting  at  the  quarter-chord.  The  differences,  while  a  small 
part  of  the  moment  about  the  midchord,  represent  the  total  moment 
about  the  quarter-chord  (aerodynamic  center),  which  is  due  entirely  to 
the  nonuniform  part  of  the  chordwise  "impressed”  velocity  distribution 
and  the  unsteady  effects;  this  is  obvious  from  equations  (43)  and  (19). 
These  moments  can,  therefore,  represent  a  significant  part  of  the 
moments  about  the  pitch  and  elastic  axes.  It  should  be  remembered 
that  the  shed  vortex  wake  is  not  well  represented  in  the  present  wake 
configuration,  and  this  could  significantly  influence  these  moments 
because  they  depend  on  higher-ordered  variations  of  the  chordwise 
"impressed"  velocity  distribution. 

Sensitivity  to  Wake  Parameters 


Computations  were  made  to  determine  the  relative  sensitivity  of  the 
results  to  variations  in  several  of  the  wake  configuration  parameters. 
The  following  were  investigated: 

1.  The  effect  of  using  a  rolled-up  root  and  tip  trailing 
vortex  versus  only  a  tip  vortex  beyond  the  grid 
portion  of  the  wake. 

2.  The  effect  of  changing  "wake  advance"  (i.  e.  ,  the 
amount  by  which  the  entire  wake  is  advanced  with 
respect  to  the  blade  producing  it). 

3.  The  effect  of  adding  distortions  of  the  tip  vortex  in 
the  Z -direction  which  approximate  that  observed  in 
smoke  studies. 

4.  The  effect  of  varying  the  radial  location  of  the  rolled- 
up  tip  vortex  (i.  e. ,  varying  the  contraction). 


For  all  computations  of  this  investigation,  the  first  three  azimuthal 
segments  of  wake  behind  each  blade  are  of  the  grid  representation 
and  the  remainder  of  the  wake  is  represented  by  concentrated  root  and/ 
or  tip  trailing  vortices. 

For  the  HU-1A  rotor  at  =  0.  26,  two  computations  were  made  to 
determine  the  effect  of  using  a  concentrated  trailing  root  vortex  at 
r/R  =  0.  30  in  addition  to  the  trailing  tip  vortex  at  r/R  =  0.  90.  The 
strengths  of  these  concentrated  trailing  vortices  were  equal  to  and 
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varied  azimushally  as  the  maximum  value  of  the  radial  distribution  of 
the  blade  bound  vcrticity.  It  was  found  that  adding  the  concentrated 
root  vortex  of  strength  equal  to  the  tip  vortex  had  a  relatively  small 
effect  on  the  airloading  and  the  bending  moments.  This  is  evident  in 
Figure  16  where  the  time  histories  of  the  airloading  at  rj R  =  0.95  and 
the  bending  moments  at  r/R  =  0.  21  are  presented.  Thus,  for  this 
rotor  and  flight  condition,  the  root  vortex  had  very  little  effect  on  the 
computed  results. 

The  effect  of  wake  advance  was  investigated  using  \.ie  HU-1A  rotor  at 
u.  -  0.  26  where  the  wake  advance  was  changed  from  0.  0  to  0.  7.  This 
had  a  relatively  small  effect,  but,  in  general,  increasing  the  wake 
advance  attenuated  and  shifted  the  phase  of  the  airloads,  blade  response, 
and  bending  moments  at  the  higher  harmonics.  This  attenuation  was 
relatively  large  for  the  response  in  the  blade  modes  near  resonance. 

For  example,  the  first  antisymmetric  blade  response  at  the  third 
harmonic  decreased  30  percent,  and  the  second  antisymmetric  response 
at  the  fifth  harmonic  decreased  38  percent.  The  attenuation  in  the 
higher  harmonics  of  the  lift  variation  is  as  expected  because  increasing 
the  wake  advance  increases  the  effectiveness  of  the  shed  vorticity 
which  attenuates  and  shifts  the  phase  of  the  lift  variations. 

The  effect  of  distorting  the  wake  in  the  Z  -direction  was  investigated  by 
introducing  distortions,  of  the  form  reported  in  references  2  and  7, 
into  the  wake  of  the  HU-1A  rotor  at  =  0.  08.  This  distortion  of  the 
tip  vortex  is  the  initial  phase  of  the  rolling  up  of  these  (initially  helical) 
vortices  into  two  wing-like  trailing  vortices  and  is  primarily  influenced 
by  the  advance  ratio,  the  disk  loading,  and  the  number  of  blades  in  the 
rotor.  The  rotor  used  to  visualize  these  distortions  and  the  HU-1  A 
rotor  are  both  two-bladed.  The  form  of  the  distortion  used  is  taken, 
therefore,  from  results  at  comparable  advance  ratios  in  reference  7 
and  scaled  according  to  the  relative  disk  loadings  for  the  two  rotors. 

It  was  hoped  that  this  would  introduce  distortions  of  at  least  the  proper 
form  and  indicate  their  effect  relative  to  the  regular  skewed  helical 
wake. 

The  distortion  of  the  wake  used  is  shown  in  Figure  17,  where  the  tip 
vortices  on  the  retreating  side  of  the  rotor  disk  are  presented.  The 
computed  airloads  at  r/R  =  0,  95  obtained  with  and  without  the  above 
wake  distortion  are  presented  in  Figure  18.  The  largest  effect  on  the 
airload  appears  on  the  retreating  side  of  the  disk  where  the  distortion 
was  large  and  such  as  to  bring  the  tip  vortex  of  the  preceding  blade 
closer  to  the  following  blade. 

The  computations  with  and  without  distortion  of  the  tip  vortices  were 
made  for  the  trailing  tip  vortices  positioned  at  r/R  =  0.  90.  In  Figure 
19  is  presented  a  comparison  of  the  airloadings  at  radial  station 
r/R  =  0.  95  when  the  same  wake  distortion  was  used  but  the  position 
of  the  trailing  tip  vortex  was  changed  from  r/R  =  0.  90  to  r/R  =  1. 00. 
The  difference  in  the  airloading  is  seen  to  be  quite  large,  and  the  results 
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are  similar  at  the  other  blade  radial  stations.  The  relatively  large 
induced  effect  due  to  the  blade  station  crossing  tip  vortices  of  the  pre¬ 
ceding  blades  (on  the  advancing  and  retreating  of  the  disk)  when  the 
vortex  is  at  r/R  -  0.  90  has  been  noticeably  diminished  when  the  vortex 
is  positioned  at  r/R  -  1.00.  This  is  also  evident  in  Figure  20  where 
the  wake-induced  velocity  variations  for  the  two  cases  are  compared. 
Similar  results  were  also  obtained  when  the  radial  position  of  the  tip 
vortices  was  changed  without  including  the  distortion  of  the  tip  vortices. 

The  magnitude  of  the  jump  in  the  wake -induced  velocity  from  the  rear 
half  of  the  disk  to  the  front  half  is  less  and  more  gradual  with  the  tip 
vortex  positioned  at  r/R  -  1. 00.  It  is  also  noted  that  the  downwash 
over  the  rear  half  of  the  disk  is  less  and  the  relative  upwash  over  the 
front  half  of  the  disk  is  greater  with  the  vortex  at  r/R  -  1. 00.  The 
reason  for  these  relatively  large  differences  in  the  induced  velocities 
and  airloads  between  the  two  cases  can  be  explained  with  the  aid  of 
Figure  21,  v/here  there  is  presented  (for  the  advancing  side  of  the  disk) 
a  projection  of  the  X-Y  plane  of  the  locus  of  blade  radial  station 
r/R  -  0.  95  together  with  the  X-Y  plane  projection  of  the  two  "youngest" 
tip  vortices.  Figur  '  21a  is  for  the  tip  vortices  positioned  at  r/R  =  0.  90 
and  Figure  21b  is  for  the  tip  vortices  positioned  at  r/R  -  1.00.  In 
Figure  21a  (for  the  tip  vortices  at  r/R  -  0.  90),  it  is  seen  that  this  blade 
station  passes  over  these  tip  vortices  at  approximately  the  intersection 
of  their  projections.  Thus  this  blade  station  passes  from  the  downwash 
side  of  the  two  vortices  to  the  upwash  side  of  both  of  them  simultaneously 
with  a  subsequent  rapid  change  in  the  induced  velocity  from  downwash 
to  upwash.*  However,  in  Figure  21b  it  is  seen  that,  for  the  vortices 
at  r/R  -  1,00,  this  does  not  happen,  but  rather  the  blade  station  passes 
over  the  two  vortices  in  a  region  where  their  projections  are  separated. 
Thus,  for  this  case,  the  blade  station  requires  a  much  larger  azimuth 
interval  to  pass  from  the  downwash  side  to  the  upwash  side'  of  both 
vortices  with  a  subsequently  slower  change  in  the  induced  velocity 
experienced.  Also,  when  the  vortices  are  at  r/R  -  1.00  (Figure  21b), 
it  is  seen  that  on  the  rear  half  of  the  disk  the  blade  station  is  farther 
from  these  two  vortices  and  on  this  front  half  of  the  disk  it  is  closer 
to  them.  Thus,  the  downwash  effect  will  be  reduced  while  the  upwash 
effect  will  be  increased  as  is  observed  in  the  induced  velocities  of 
Figure  20. 


It  should  be  remembered  when  viewing  these  wake-induced  velocity 
time  histories  with  respect  to  the  effects  of  the  preceding  blade  wakes, 
that  they  represent  the  total  wake-induced  velocity  experienced  by 
the  blade  station  and  include  the  downwash  of  the  attached  near  wake 
trailing  vorticity  which  superimposes  a  negative  bias  on  the  effects 
of  the  preceding  blade  wakes. 


Sensitivity  to  Blade  Natural  Frequencies 


For  the  HU-1A  at  ll  =  0.  26  the  computed  response  of  the  second  anti¬ 
symmetric  blade  bending  mode  was  overpredicted  (see  Figures  8c  and 
8d),  and  this  mode  is  very  near  resonance  with  the  fifth  harmonic  of  the 
rotor  speed.  This  situation  occurs  for  most  of  the  computed  cases 
because  it  is  impossible  to  avoid  proximity  to  resonance  for  all  the 
blade  modes.  Thus,  because  the  damping  ratios  are  quite  small  (it  is 
estimated  to  be  0.  035  from  the  computed  response  and  phase  angle  of 
the  second  antisymmetric  mode  of  the  HU-1A  at  /j.  -  0.  26),  the  operating 
point  for  these  bending  modes  can  be  on  the  very  steep  slopes  of  the 
amplification  factor  curve  and  small  changes  (or  errors)  in  the  natural 
frequency  can  have  very  large  effects  on  the  response  in  such  modes. 
This  was  confirmed  by  making  small  changes  in  the  natural  frequency  of 
the  second  antisymmetric  mode  for  the  HU-1  A  at  /c  =  0.  26.  It  should 
be  noted  that  this  is  not  a  characteristic  of  the  computational  procedure, 
but  rather  an  inherent  characteristic  of  the  physical  system  and  will 
be  a  problem  for  any  method  of  prediction.  Because  of  this  sensitivity, 
the  blade's  coupled  natural  frequencies  should  be  used  even  if  the 
uncoupled  mode  shapes  are  used. 


CONCLUSIONS  AND  RECOMMENDATIONS 


The  results  of  this  effort  have  proved  the  feasibility  of  this  computational 
approach  to  the  rotor  aeroelastic  problem.  The  following  objectives 
have  been  attained: 

1.  "Closing  the  loop1'  with  respect  to  the  blade  motions. 

2.  Satisfying  the  chordwise  aerodynamic  boundary  condition. 

3.  Simplifying  the  discrete  wake  representation. 

Considering  the  sensitivity  of  the  results  to  the  blade  natural  frequencies 
and  the  wake  configuration,  the  agreement  of  the  computed  results  with 
the  measured  results  is  quite  good. 

Based  on  this  work,  the  following  specific  conclusions  are  made 

1.  The  computed  bending  moments  and  airloads  are 
sensitive  to  the  blade  natural  frequencies. 

2.  Accurate  computation  of  the  blade  airloads  and  bending 
moments  requires  the  coupled  natural  frequencies  of  the 
blades  to  be  accurately  determined. 

3.  The  computed  results  are  sensitive  to  the  wake  configur¬ 
ation.  For  the  discrete  type  wake  representation  used, 
it  was  found  that,  for  the  cases  computed, 

a.  The  presence  of  a  concentrated  root  vortex  had  little 
effect. 

b.  Advancing  the  shed  vortices  (closer  to  the  blade 
trailing  edge)  attenuated  and  shifted  the  phase  of 
the  computed  blade  response. 

c.  Distortions  of  the  wake  normal  to  the  tip-path 
plane  can  have  a  significant  influence  on  the  results. 

d.  The  position  of  the  rolled-up  tip  vortex  can  have  a 
significant  effect. 

4.  Significant  aerodynamic  pitching  moments  about  the 
quarter-chord  are  computed  (for  a  flat  plate  or 
symmetrical  airfoil)  by  satisfying  the  'chordwise  aero¬ 
dynamic  boundary  condition  and  including  unsteady  effects. 


The  present  procedure  does  not  impose  the  rotor -fuselage  trim  con¬ 
straints  on  the  resultant  rotor  forces  and  moments;  that  is,  the  rotor 
shaft  angle  and  rotor  control  settings  are  required  as  inputs  to  the 
computation  and  the  rotor  resultant  force  and  moment  are  unconstrained 
results  of  the  computation.  However,  in  the  task  of  designing  a  rotor 
system  for  a  specific  fuselage,  c.  g.  location,  and  flight  condition,  it 
is  the  rotor  control  settings,  blade  stresses,  and  aerodynamic  loads 
which  are  sought  and  the  resultant  rotor  forces  and  moments  must 
balance  the  forces  and  moments  of  the  fuselage.  Thus  for  the  design 
and  development  of  rotor  systems,  and  also  for  the  analysis  of  existing 
rotor  systems,  a  practical  method  of  prediction  should  include  the 
appropriate  trim  equations. 

The  computational  procedure  in  its  present  form  is  incomplete  and 
has  limited  use  for  rotor  design  analysis.  On  the  basis  of  the  encour¬ 
aging  results  which  have  been  obtained  with  the  present  method,  it  is 
recommended  that  further  effort  be  expended  to  remove  its  limitations. 
Specifically  the  following  should  be  accomplished; 

1.  Incorporate  the  pitch,  torsion,  and  in-plane  degrees  of 
freedom  and  the  significant  coupling  between  all  degrees 
of  freedom. 

2.  Develop  a  procedure  for  adequately  predicting  the 
wake  configuration. 

3.  Incorporate  the  necessary  rotor -fuselage  trim  equations 
into  the  procedure. 
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BLADE 


Figure  2.  THE  TIP-PATH-PLANE  ORIENTED  COORDINATE  SYSTEM  USED 
FOR  COMPUTING  THE  INDUCED  VELOCITIES. 


Figure  3.  EXAMPLE  OF  THE  K-SUBSCRIPT  NOTATION  FOR  AIRLOAD 
POINTS  IN  ROTOR  DISK  WHEN  NR  =  3,  NA  =  8. 
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Figure  6.  MAJOR  STEPS  IN  THE  ITERATIVE  PROCEDURE  OF  SOLVING  THE  EQUATIONS 
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a.  TIME  HISTORIES  OF  LIFT 

Figure  8.  MEASURED  AND  COMPUTED  LIFT  AND  BENDING  MOMENTS 
HU- 1 A  AT/x  =  0.26. 


b,  HARMONICS  OF  LIFT 


Figure  8  (Cont'd).  MEASURED  AND  COMPUTED  LIFT  AND  BENDING  MOMENTS 

HU- 1 A  AT  fx  =  0.26. 
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c.  TIME  HISTORIES  OF  BENDING  MOMEHTS 


Figure  8  (Cont'd).  MEASURED  AND  COMPUTED  LIFT  AND  BENDING  MOMENTS  ; 

HU- 1 A  AT /z  =  0.26. 
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Figure  9 


a.  TIME  HISTORIES  OF  LIFT 

MEASURED  AND  COMPUTED  LIFT  AND  BENDING  MOMENTS 
Hll-  I A  AT  =  0,08. 
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d.  HARMONICS  Of  BENDING  MOMENTS 


Figure  9  (Cont'd).  MEASURED  AND  COMPUTED  LIFT  AND  BENDING  MOMENTS; 

HU-  i  A  AT//.  =  0.08. 
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Figure  10a.  COMPUTED  OUAS i -STEADY  rt  TOTAL  r,  AND 
LIFT  AT  r/R  =  0.95;  HU-IA  AT/*=  0.08 
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b.  HARMONICS  OF  LIFT 

Figure  II  (Cont’d).  MEASURED  AND  COMPUTED  LIFT  AND  BENDING  MOMENTS 

K-34  AT*.  =  0.29. 
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c.  TIME  HISTORIES  OF  BENDING  MOMENTS 

Figure  I!  (Cont'd)  MEASURED  AND  COMPUTED  Li  FT  AND  BENDING  MOMENTS; 

H-34  AT/a=  0.29. 
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a.  TIME  HISTORIES  OF  LIFT 

Figure  12.  MEASURED  AND  COMPUTED  LIFT  AND  BENDING  MOMENTS 
H~34  AT/*  =  0. 18. 
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c.  TIME  HISTORIES  OF  BENDING  MOMENTS 

Figure  \2  (Cont'd)  MEASURED  AND  COMPUTED  Li  FT  AND  BENDING  MOMENTS; 

H-34  AT /a  =  0. 18. 
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Figure  14a.  COMPUTED  PITCHING  MOMENTS;  HU-IA  AT  /ju  =  0.26, 
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Figure  14b.  COMPUTED 
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Figure  15,  COMPUTED  PITCHING  MOMENTS  ABOUT  MID-CHORD;  ti-M  AT  pc  -  0.18. 
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gure  13.  COMPARISON  OF  LIFT  AT  r/R  =0.95  COMPUTED  WITH  AND 
WITHOUT  WAKE  DISTORTION.  HU-IA  AT  jjl  =  0.08. 
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Figure  19.  COMPARISON  OF  LIFT  AT  r/R  =  0.95  COMPUTED  FOR  TWO  POSITIONS 
OF  THE  TIP  VORTEX  IN  THE  WAKE.  HU-IA  AT  jjl  =  0.08. 
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20.  COMPARISON  OF  THE  INDUCED  VELOCITY  AT  r/R  =  0.95 
FOR  TWO  RADIAL  POSITIONS  OF  THE  TIP  VORTEX  IN  THE 
WAKE;  HU- 1 A  AT  ,u  =  0.08. 
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Figure  21.  RELATIVE  ORIENTATION  OF  THE  BLADE  AND  RECENT  TIP  VORTICES 
ON  ADVANCING  SIDE  OF  THE  DISK;  HU- I A  AT  =  0.08. 


APPENDIX  I 


TWO-DIMENSIONAL  OSCILLATING  WING  INVESTIGATION 
OF  WAKE  SHED  VORTEX  REPRESENTATION 


A  numerical  computational  procedure  was  developed  and  programmed 
for  computing  the  Lift  and  pitcning  moments  experienced  by  a  two- 
dimensional  oscillating  wing  whose  shed  wake  was  represented  by 
equally  spaced  concentrated  vortices.  This  procedure  is  parallel  to 
that  developed  and  used  to  represent  the  rotor  blade  and  described 
under  the  section  Blade  Representation.  For  this  problem,  however, 
the  integrals  for  the  cosine  components  of  chordwise  induced  velocity 
distribution  due  each  wake  vortex  element  can  be  evaluated  in  closed 
form.  The  procedure  is  outlined  below,  and  the  results  of  the  investi¬ 
gation  are  presented  and  discussed. 

The  notation  to  be  used  is  illustrated  in  the  following  diagram 


The  airfoil  chord  is  centered  on  the  x  -axis  with  the  leading  and 
trailing  edges  at  x  =  -1  and  Y  -  +1  respectively.  The  airfoil  shed 
vorticity  is  represented  by  a  distribution  of  concentrated  equally 
spaced  two-dimensional  vortices  and  the  nondimensional  spacing 
given  by 

7  ZttU 

d.  -  — (half-chords  traveled  per  shed  vortex)  (I _  1 ) 

X  b  <j* 

where 

U  =  freestream  velocity 

b  -  airfoil  semichord 

u)  -  frequency  of  airfoil  oscillations 
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X 


=  number  of  shed  vortices  per  cycle  of  the  airfoil  motion 
used  in  the  wake  representation 


The  nondimensional  distance,  O  t  from_the  midchord  to  the  first  shed 
vortex  was  taken  to  be  independent  of  d  .  When  D  =1  +  d  ,  the  first 
shed  vortex  is  behind  the  trailing  edge  a  distance  equal  to  the  wake 
spacing  and  conceptually  the  vortices  are  being  released  from  the 
trailing  edge  for  the  airfoil  in  each  position,  h  ,  of  its  oscillation. 
When  D<1  +  d  >  the  wake  is  considered  to  be  advanced  with  respect  to 
the  airfoil.  The  distance  of  each  shed  vortex  from  the  midchord  is 


Sn  ~  (n-f)  d  +  D  (1-2) 

and  n  is  the  relative  position  number  of  the  vortex  behind  the  trailing 
edge.  The  airfoil  position  number,  k,  ,  is  a  label  to  indicate  the 
relative  position  of  the  airfoil  in  its  cycle  of  oscillation  (  fc  =  1  is  the 
beginning  of  the  oscillation  and  there  are  X  positions).  The  values  of 
the  airfoil  bound  vorticity,  rk  ,  at  each  position  in  the  cycle  of 
oscillation  are  the  unknowns  and  because  of  the  periodicity  of  the 
problem  they  repeat  at  corresponding  position  numbers;  there  are  thus 
X  unknown  values  of  /*  .  The  strengths  of  the  shed  vortices  and  their 
positive  sense  are  indicated  in  the  diagram  for  the  first  few  vortices 
behind  the  airfoil  as  the  difference  in  strength  of  the  bound  vorticity 
between  successive  airfoil  positions. 

The  chordwise  aerodynamic  boundary  condition  on  the  airfoil  is 
expressed  by  the  following  integral  equation, 


(  -  h  +  Ucc  +  a  ba) 

n=t 


*rn 

2ir(Zn-x)b 


/  r+f  y(x)  d.% 

2tt  i /  (x-x) 


d-3) 


where  A/'  is  the  number  of  cycles  of  v/ake  representation  retained  and 
h  ,  cx  ,  and  cx.  are  the  blade  input  oscillation  in  plunge  and  pitch.  By 
letting  x  =  -cose  and  expanding  the  bounc  vorticity  distribution  in  a 
Glauert  series 


yC&)  -  2 (Aa  coi  j  +  22  An  sin  nG  )  , 
*  n=f 


d-4) 


the  integral  of  (1-3)  can  be  evaluated  as  the  following  cosine  series: 


J+1  Y(x}dx_ 


OO 


Zrr  (x-x) 


-  ~  A0  +  JT  A„  cos  . 


d-5) 


n*t 
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Thus  equation  (1-3)  can  be  written  as 


oo  "Z.  Ar„/b 

-h+U*-bacose-A0  +  Z.A„cosne+Z  Zjr(^~co,e) 


0.  (1-6) 


By  expanding  each  term  of  the  summation  in  a  cosine  series  and  then 
equating  coefficients  of  like  terms  in  the  resulting  equation,  the 
following  expressions  for  the  Glauert  coefficients  are  obtained: 


7TbAa = 

-^2irb(-h  +■ 

NT 

zZ  *r„N0„ 

*  ns  i 

NX 

rrbA1  * 

-rrb2<k  -£ 

n-t 

VX 

TrbA^  - 

“Z  &Cn  ;  ^ 

=  234  •  •  ■ 

9  9  *  % 

a*/ 


where 


(1-7) 


urn-f  -*n) 


m 


(1-8) 


is  the  result  of  evaluating  the  Fourier  cosine  coefficients  of  the  terms 
of  the  summation  of  (1-6). 


The  total  bound  vorticity  for  blade  position  k  is  obtained  by  integrating 
y(@)  (equation  1-4)  over  the  chord.  Thus, 

r k.  s  2irb  (Aa  -h  j A/)^ 


where 


AX 


n  =  it  +Z 

n=t  ' 

NT 

rk  =  ik  +Z  Arn  Ckn 


Ik  =  \zrrb  (-h  +  Uoc,)  +  vb2  6c, 


d-9) 


is  the  quasi-steady  part  of  rk  ,  and 

£kn  B  (^on  ’  n  )fc  • 

Because  the  Afn  are  differences  of  the  X  values  of  PK  ,  (1-9)  can  be 
written  as 
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(1-10) 


rk  ~  Ik  *22  ~^k)l  Xk.1  *  k-Z  ~  T c-t)z  Ck2  *  '  "  }• 

n-f  i  J 

By  regrouping  the  summation  on  like  ,  a  set  of  X  simultaneous 
equations  in  are  obtained  as 


X 


/I  *  *  Z  <ny  0 

j 


(i-ii) 


where  <7£  is  the  summation  of  all  which  in  (I - 1 0)  are  coefficients 

of  r-  .  J 


For  specified  I*  ,  the  set  of  equations  (I - 1 1 )  are  solved  for  the  rK  and 
then  the  Glauert  coefficients  evaluated  by  equations  (1-7).  With  the 
Glauert  coefficients,  the  ratios  of  lift  and  moment  to  their  quasi-steady 
values  are,  for  each  position,  k  ,  of  the  airfoil,  computed  fron. 


Z7rbk 

Ik 

2rrbk 

Ik 


Ac  + 
Aa  + 


k 

k 


d-12) 


where 


L  o  -  /?  UIk 
M0  =  />UIk*±. 


(1-13) 


The  Glauert  coefficients  are  complex  and  thus  the  L/.  and  /M  will  be 

*  'uO  ,nO 

complex. 


The  relative  merit  of  this  wake  representation  was  investigated  by 
computing  the  lift  and  pitching  moment  transfer  functions  (i.  e.  ,  L/L0 
and  M//*fc  versus  reduced  frequency,  bu>/U  )  and  comparing  with  the 
classical  results  as  presented  by  Karman  and  Sears  in  reference  8  for 
an  airfoil  oscillating  in  plunging  or  in  pitching.  It  should  be  noted  that 
the  results  as  presented  for  oscillation  in  pitch  represent  only  the  pitch 
rate  effects  (i.  e.  ,  due  to  a.  only);  this  was  done  because  the  effect  of 
the  pitching  displacement  is  the  same  as  for  plunging  velocity  (i.  e.  , 
both  h  and  <x  contribute  only  to  the  constant  component  of  the  impressed 
chordwise  velocity  distribution  and  thus  affect  only  A0  whereas  oc  affects 
only  Af  . 

The  quasi-steady  component,  Ik  ,  of  the  bound  vorticity  depends  on  both 
the  pitching  and  plunging  motions  of  the  blade.  Thus  for  the  computations 
(whether  for  pitching  or  plunging  oscillations),  Ik  is  always  set  equal 
to  !-c.osu)t  and  then  interpreted  as  either  the  quasi-steady  part  of  A0  or 
Af  in  evaluating  the  Glauert  coefficients.  That  is,  the  difference  in  the 
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computations  for  the  plunging  and  pitching  cases  is  only  in  the  evaluation 
of  the  A0  and  A,  Glauert  coefficients.  Because  Ik  was  taken  as 
i-cosaJi  ,  the  cosine  components  of  the  resulting  L/L0  and  M/M0  are 
interpreted  as  the  real  parts  and  the  sine  components  as  the  negative 
imaginary  parts. 

In  Figure  I.  1  are  presented  the  results  of  the  computations  ’hen  the 
wake  advance  was  set  equal  to  zero  (i.  e.  ,  D  -  /  +  ).  That  is,  the 

distance  of  the  first  shed  vortex  behind  the  trailing  edge  is  always  equal 
to  the  shed  vortex  spacing,  oL  .  The  computations  were  made  for  X  = 

3.  8,  30,  and  90  shed  vortices  per  cycle  of  the  blade  oscillation  and  the 
reduced  frequency,  buJ/U  range  from  0  to  0.  8.  Figure  I.  la  presents  the 
L/L0  and  for  the  airfoil  oscillating  in  plunge,  and  Figure  I.  lb 

presents  them  for  the  airfoil  oscillating  in  pitch.  The  agreement  with 
t.he  classical  solution  improves  as  X  increases,  which  is  to  be  expected 
because  the  wake  representation  is  approaching  the  continuous  wake.  A 
vector  drawn  from  the  origin  to  a  point  on  one  of  the  curves  represents 
the  magnitude  and  phase  angle  {relative  to  the  quasi -steady  values)  of 
that  quantity  for  reduced  frequency  corresponding  to  the  point;  note  that 
the  scale  of  the  abscissa  has  been  broken  between  0  and  0.  5  on  several 
figures.  It  is  noted  that  the  differences  between  the  approximate  wake 
results  and  the  exact  wake  results  is  primarily  in  the  phase  and  that 
the  agreement  is  better  for  the  pitching  moments  than  the  lift. 

The  computations  were  repeated  for  several  values  of  the  wake  advance 
relative  to  the  blade  and,  although  the  results  improved  with  increasing 
wake  advance,  the  optimum  appeared  to  be  about  70  percent.  The 
results  for  the  wake  advanced  70  percent  (i.  e.  ,  D  -  1  +  d/3  )  are  pre¬ 

sented  in  Figure  I.  2.  It  can  be  seen  that  the  curves  for  various  numbers 
X,  of  shed  vortices  per  cycle  are  grouped  much  closer  to  ^he  classical 
results  and  that  the  discrepancies  in  both  the  magnitude  and  phase  angle 
have  been  reduced.  The  principal  difference  between  the  approximate 
and  classical  results  is  still  in  the  phase  angle;  the  pitching  moments, 

M/ Mg  ,  are  in  better  agreement  than  the  lift  L/L0 . 

It  should  be  remembered  that  the  curves  presented  from  this  investiga¬ 
tion  are  for  a  constant  number,  X  ,  of  shed  vortices  per  cycle.  The 
use  of  a  similar  wake  representation  for  a  rotor  does  not  result  in  X 
being  a  constant,  but  rather  X.  is  then  a  function  of  the  blade  radius 
and  harmonic  order  of  the  blade  dynamic  response.  For  the  rotor 
wake  application,  it  is  the  shed  vortex  spacing,  <£  ,  which  is  relatively 
constant  (there  is  some  variation  with  azimuth  angle,  *jJ  ,  from  the 
retreating  to  the  advancing  sides  of  the  rotor)  and  the  number,  X  ,  of 
shed  vortices  per  cycle  decreases  with  increasing  harmonic  order 
(i.  e.  ,  the  lowest  harmonics  have  the  best  representation). 
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0  0.5  0.6  0.7  0.8  0.9  1.0 

LIFT,  l/L0  -  REAL  PART 


PITCHING  MOMENTS,  M/M0  -  REAL  PART 
a.  PLUNGING  OSCILLATIONS 


Figure  I. I  TWO-DIMENSIONAL  OSCILLATING  AIRFOIL  LIFT  AND  PITCHING  MOMENT 
TRANSFER  FUNCTIONS;  ZERO  WAKE  ADVANCE  (D  =  l+d). 


77 


PITCHING  MOMENTS, 

M  -  IMAGINARY  PART  LIFT,  L/L_  -  IMAGINARY  PART 


s.  PLUNGING  OSCILLATIONS 

Figure  1.2  .  TWO-DIMENSIONAL  OSCILLATING  AIRFOIL  LIFT  AND  PITCHING 
MOMENT  TRANSFER  FUNCTIONS;  0.70  WAKE  ADVANCE  (D  =  |+d 
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PITCHING  MOMENTS, 

L  -  IMAGINARY  PART 


APPENDIX  II 


THE  CHORDWISE  INDUCED  VELOCITY  DISTRIBUTION 
DUE  TO  A  WAKE  VORTEX  ELEMENT 


The  notation  to  be  used  is  defined  in  the  following  diagram,  where  CL  and 
b  are  the  end  points  of  the  straight-line  vortex  filament  and  p  is  the 
point  at  which  the  induced  velocity  is  to  be  computed.  These  three 
points  can  be  arbitrarily  oriented  in  space. 


P 

The  differential  induced  velocity  at  a  point  due  to  an  elemental  length 
of  vortex  filament  is  given  in  vector  notation  as 

<u 

This  can  be  put  in  the  following  form: 


(II-2) 
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The  cosines  of  the  angles  oc.  ,  /}  ,  and  /  subscripted  with  a>  and  5  are, 
respectively,  the  direction  cosines  of  the  vectors  F  and  of  the 
diagram  on  the  previous  page. 

The  total  velocity  at  p  due  to  the  vortex  filament  is  obtained  by  integra¬ 
ting  (11-2)  over  the  length  of  the  filament  from  a  to  b  \  the  result  is 


r h  r  cL  -  £/a. -  rj, )  ecs  &g, 
?  ~  *7T  L rlrbCl-cos'&cJ) 


K  . 


(II-4) 


For  use  in  the  numerical  computations,  this  expression  for  the  velocity 
is  written  in  terms  of  the  direction  numbers  (and  thereby  the  coordinates 
of  the  three  points)  of  the  two  radius  vectors  from  the  filament  end 
points  to  the  point,  p  ,  as  follows: 


(  r/  +tt )  ( ro- *rb) 


rh  +  +  n*  +  n^nb) 


(  (^nb  +  (n_5) 


(-4.  nb  +J,b  +  (  £a,  mh  rr?^)  k  J 


where 


t 

ni.  =  f  4 


w 


/V  -  A  *  +  ^  y 


>  L  =  a>  or  b 


are  the  direction  numbers  and  magnitudes  of  the  vectors,  ^  an<^  r!>  • 
The  (z.y,  $  )  are  the  coordinates  of  the  point,  p  ,  and  (  y  >  £  )  ar:  the 
coordinates  of  the  vortex  filament  end  points  with  the  corresponding 
subscript. 


Only  the  £  -component  of  the  chordwise  induced  velocity  distribution  is 
used  in  the  present  method.  The  expression  for  this  chordwise  distri¬ 
bution  can  be  derived  from  the  |  -  component  of  (II - 5)  by  considering  it 
to  express  the  induced  velocity  at  any  one  point  (  x,.  y ,  ^)  on  the  chord 
and  expressing  this  point  in  terms  of  a  chordwise  coordinate  S  (relative 
to  the  midchord)  by  the  following  expressions: 

Ax 

X  =  -  S’ 

m  b 


y  ’ 


t  = 


(H-6) 
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where 


Ax  =  xL  -  xM 
Atj  -  Vt.  ~  V* 

A  $  *  H  ~ 


are  the  direction  numbers  for  the  semichord  (subscripts  L  and  M 
designate  the  leading  edge  and  the  midchord  points,  respectively).  When 
these  expressions  are  substituted  for  x  and  j  in  the  j  -components 
of  (II-5),  the  following  expression  for  the  -component  of  the  chord- 
wise  induced  velocity  is  obtained  in  terms  of  the  chordwise  coordinate  s 


(w)  yc,  +  c2s  +  c3s*  +Css  +  Cti\c7  +  C0&  +  C9s*] 


^  ijc,  +cts+c3s*  +  ->c3s2  -jc^*CsS  +  Cfcs*  +Ci0+C,,s  +  Cns* 


(II -7) 


where 

C1  -  (xH  -gtJ}  +  (<Jtf  ~  fa)  ''  1  ~  fa.) 

-€*.)&*+(!/* -*1*)*#  + 
c3  =  p  | A%*  +  Ayz  +  A^  =  / 
c*  =  (x-rt  "ft,)  +  +  ~  £i) 

cs  =  "  f  4*  *(}■**- &) 

C(,  -  -  / 

(II  -  8) 

C7  =  (.**  -£*,)(</„  -yb )  - 

Cg  -  -  J  [{C**“£j  ~C*/V  {(VM-ft,  }  - 

c9  -  o 

Cio  *(X*-$c*)(Xrt'4i)  +  (*M'-'?a,)(y*->?b)  +(i-M -cJ(?*  -?*) 

Cff  =  £  (cz  +  Cf  J 

CfZ  -  ^3  -  t  • 
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For  the  present  application,  the  first  four  Fourier  cosine  components 
of  the  induced  velocity  distribution  given  by  (II- 7)  are  desired  but  the 
integrals  for  the  Fourier  coefficients  of  (11-7)  could  not  be  evaluated. 
The  complexity  of  equation  (II-7)  arises  because  it  expresses  the  chord- 
wise  induced  velocity  fc~  all  relative  orientations  of  the  chord  and 
vortex  filament,  including  the  situation  where  they  are  extremely  close 
(or  intersecting)  and  the  cnordwise  velocity  variations  are  very  rapid. 

It  was  assumed  that,  for  this  application,  the  chord  to  wake  filament 
distances  encountered  in  the  computer  program  will  never  be  extremely 
small  and,  therefore,  the  corresponding  chordwise  variations  of  the 
induced  velocity  due  to  the  elements  will  not  be  very  rapid.  Thus  when 
(II  — 7)  is  expanded  in  a  power  series,  about  s  =  0,  only  the  first  four 
components  are  retained  and  this  finite  power  series  is  transformed  to 
a  cosine  series.  The  resultant  power  series  is: 


4-7T< 


?  - 


\ot,  oc.0Ai 

_A,  Ao*  J 

l  ^ 

+ 

Z  GZj&f 

*  a>0(£f  -  **  Y 

.  Ao  Ac  2 

VA  Ao  I. 

(H-9) 
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Ao  ft*  'U2 


A2 


V  -  .4l  -  43)1 

)  (a3  Ac  PoK 


where 


8  ^7  (Vq  +  -f c If  j 

=  c’[fc(ic 

**  *  i-i%) MtS •?  M 

(it  -  W } 


oe. 


OC,y  = 


<ZC,  8 


fF) 


+  i  £ 


Po a  c  f  *  c  f0  y Af 

A,  =  CtCr  +  Cz  c*  +  T&r«[c„  +  ” ) 
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/ - r  /  /  i  c£  i  czcs  /  /  c/  '| 

A  =  c,*c,c^r^y^[oo(— +  t-F^+YFy-77T) 

+  £il(S2-  +  £l)  +  fl 

2  \Cj  C*  )  J 
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(11-10) 


and  the  C‘s  are  defined  in  (II-8).  In  the  computer  program  a  further 
hierarchy  of  symbolization  was  used  to  take  advantage  of  the  repetition 
of  certain  groups  and  combinations  of  the  C's  and  thereby  increase 
the  computational  efficiency. 

The  transformation  from  the  power  series  (II -9)  to  the  following  cosine 
series 


H-7T 

-JT  ~  S0  +  5,  COS0+  sz  COS  ZG  +s3  cos 


(Il-il) 


(11-12) 


is  effected  by  the  following  transformation: 

s  -  -  b  cos  0 

where 

Sf  3  _cv*  ~~*a'3  (11-12) 

si  -  "z  Q't 
5S  ~  ~  O' 3 

and  the  an{n  -  0,  1,  2,  3)  are  the  coefficients  of  (II-9)  multiplied  by 
b  .  Thus  (11-12)  are  the  desired  expressions  for  the  first  four  cosine 
coefficients  of  the  chordwise  induced  velocity  distribution  due  to  a 
single  wake  element,  but  it  should  be  noted  they  have  been  approximated 
by  truncating  the  power  series. 
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APPENDIX  III 


EFFICIENT  METHOD  FOR  NUMERICALLY  COMPUTING  DERIVATIVES 

OF  PERIODIC  FUNCTIONS 


If  a  periodic  function  is  numerically  defined  at  MA  equally  spaced 
points  on  its  period,  it  is  possible  to  represent  it  with  a  finite  Fourier 
series  of  NA  terms  (reference  12,  page  430).  When  NA  is  an  even 
number,  the  series  is 


where 


NA 

*(?k)=  £  cos  n  +£  bn  sin  n  *k 

O ,0  n~  f 


,  k.zNA 

---E 


*(*,) 


2  \~  I 

k  =  / 

2  k-t/4 

bo  -  -j-T  f(fk)sinn^ \ 

h* //A 

°-na  =  *(A)c*sir(k‘t) 

■y  k-t 

yk  -(k'-l)A'P 

W  2rr/„„ 


n  4  o 

n<  OC 
Z 


(III-l) 


(III -la) 


(III -lb) 


(III-lc) 


(HI- Id) 


When  NA  is  an  odd  number,  the  last  cosine  term,  a 


"A/z 


is  absent. 


The  first  derivative  of  the  function  can  also  be  represented  by  a 
Fourier  series  as  follows: 


fSA  NA  j 

*(+k)  -  -t.  na,„  sir?  n  (>k  cos  n\j/k  . 


(IH-2) 


n  =  / 


n  &  f 
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If  the  expressions  (III - 1 )  are  substituted  for  the  coefficients  of  (III -2) ,  the 
derivative  can  be  written  as 


rjA  _/ 

-Z  • k) 


(III -3) 


which  expresses  the  derivative  at  a  point,  ^  ,  in  terms  of  the  magnitude 

of  the  function  at  the  NA  points,  ^  .  It  is  noted  that  the  expression  in 

braces  is  only  a  function  of  the  number,  NA  ,  of  points  at  which  the 
function  is  defined  and  thus  can  be  precomputed  as  a  set  of  coefficients. 

This  procedure  is  equivalent  to  solving  for  the  Four  er  coefficients  of 
the  function  by  equation  (III - 1 )  and  then  evaluating  the  derivatives  by 
(III-2);  however,  it  eliminates  the  step  of  evaluating  the  Fourier  coeffi¬ 
cients  by  precomputing  coefficients  of  tffl)  in  (III - 3) .  This  will  save 
a  considerable  amount  of  numerical  computation  when  the  derivatives 
must  be  evaluated  for  more  than  one  function  and  each  function  is 
defined  at  the  same  number,  NA  ,  of  points. 

This  procedure  was  used  in  the  present  prediction  method  because  the 
derivatives  of  two  combinations  of  Glauert  coefficients  must  be 
evaluated  for  each  radial  station  (thus  2  ~  NR  functions)  and  this  must 
be  repeated  for  each  iteration. 

It  should  be  noted  that  this  procedure  is  also  applicable  to  the  higher 
derivatives.  For  example,  the  expression  for  the  2nd  derivative  is 


ha  f 

■(M  +  t 

/  c  n*f  L  J  v 


(III-4) 
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APPENDIX 


IV 


FLOW  DIAGRAM  FOR  THE  COMPUTER  PROGRAM 

(PARTS  1  AND  2) 


These  flow  diagrams  do  not  present  all  the  details  of  the  program  but 
only  the  general  flow  of  both  parts  and  their  subroutines.  The  flow 
diagrams  for  program  Parts  1  and  2  and  their  subroutines  are  pre¬ 
sented  in  Figure  IV.  1  and  IV.  2,  respectively.  The  symbols  used  in 
these  diagrams  are  defined  below. 


COORD  Name  of  subroutine  which  computes  blade  and 
wake  coordinates 

CVLM  The  convergence  limit  for  both  the  r's  and  the 
generalized  coordinates 

GHCRD  Name  of  subroutine  which  solves  the  equations  of 
motion  for  the  response  in  the  generalized 
coordinates 

GNFRCE  Name  of  subroutine  which  computes  the  lift, 
pitching  moments  and  generalized  forces 

X  Constant  term  of  the  "  r  "  equation 
{-ft  +  V,  <X9  +  -jbcz 5) 

I8LD  Blade  number 

IRTR  Rotor  'ition  number 

IS  Mode  (or  degree  of  freedom)  number 

1ST  Giauert  coefficient  subscript 

I3YM  Control  word  which  indicates  whether  the 

mode  under  consideration  is  symmetric  (=  +1) 
or  antisymmetric  (=  -1) 

UG  Iteration  number  for  the  r  -equations 

jtX  Iteration  number  for  the  blade  motions 

ITYP  Control  word  v/hich  indicates  whether  the  rotor 

has  a  teetering  hub  (=  1)  or  a  rigid  hub  (=  2) 
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IWN  Blade  number,  whose  wake  is  under  analysis 
^  Blade  spanwise  segment  number 
Azimuth  position  number 
Number  of  azimuth  positions 
Number  of  blades 


ND  Number  of  blade  modes  (degrees  of  freedom) 

HER  Number  of  equally  spaced  radial  blade  segments 

at  which  the  computed  lift  is  interpolated.  This  is 

for  use  in  evaluating  integrals  of  the  generalized  airload. 

NrTG  Maximum  number  of  iterations  allowed  on  the  P  -equations 

HUM  Maximum  number  of  iterations  allowed  on  the  blade 

motions 


HR  Number  of  blade  spanwise  segments 

RA  Distance  from  wake  vortex  element  to  blade 
segment 

RTOGAM  Ratio  of  the  change  (between  successive  iterations) 
in  r  to  its  present  value 


RTOGCO  Ratio  of  the  change  (between  successive  iterations)  in 
tip  deflection  (for  each  mode)  to  its  present  value 

ZH1 ,  ZN2,  ZH3  Zones  of  influence  about  a  blade  segment.  They 

determine  how  the  induced  effects  (due  to  each  wake 
element)  on  the  blade  segment  are  to  be  computed. 
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a.  MAIN  PROGRAM  OF  PART  I 

Figure  U.  I  FLOW  DIAGRAM  FOR  PART  I  OF  COMPUTER  PROGRAM. 


90 


c.  SUBROUTINE  COORD  OF  PART  i 


Figure  U.  I  (Cont'd)  FLOW  DIAGRAM  FOR  PART  I  OF  COMPUTER  PROGRAM. 
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a.  MAI H  PROGRAM  OF  PART  2 

Figure  IT.  2  FLOW  DIAGRAM  FOR  PART  2  OF  COMPUTER  PROGRAM. 
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INDICATE  r'x  FAILED 
TO  CONVERGE  AFTER 
NITG  ITERATIONS 


b.  SUBROUTINE  GAMMA  OF  PART  2 


Figure  IT.  2  (Cont'd)  FLOW  DIAGRAM  FOR  PART  2  OF  COMPUTER  PROGRAM 
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c.  SUBROUTINE  GNFRCE  OF  PART  2 


Figure  U.2  (Cont’d) 


FLOW  DIAGRAM  FOR  PART  2  OF  COMPUTER  PROGRAM 
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d.  SUBROUTINE  GNCRD  OF  PART  2 

gure  U.  2  (Ccnt'd)  FLOW  DIAGRAM  FOR  PART  2  OF  COMPUTER  PROGRAM. 
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